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Development of a Linearized Unsteady Euler Analysis 
for Turbomachinery Blade Rows 


Summary 


Progress toward the development of a linearized, unsteady, aerodynamic analysis for 
axial-flow turbomachinery blading is outlined in this report. The linearization is based on 
the Euler equations of fluid motion and is motivated by the need for an efficient aerodynamic 
analysis that can be used in predicting the aeroelastic and aeroacoustic responses of blade 
rows. The field equations and surface conditions required for nonlinear and linearized, invis- 
cid, unsteady, aerodynamic analyses of three-dimensional flows through single, stationary or 
rotating, blade rows operating within cylindrical ducts are derived. Also, a description of the 
approximate far-field conditions that are currently available for completing these unsteady 
aerodynamic formulations is given. 

An existing numerical algorithm for determining time-accurate solutions of the nonlinear 
unsteady flow problem is described, and a numerical model, derived from this nonlinear flow 
solver, is formulated for the first-harmonic, linear, unsteady problem. In previous work at 
NASA Lewis, two-dimensional versions of the nonlinear aerodynamic and numerical models 
have been implemented into a time-accurate, unsteady flow code, known as NPHASE. Under 
the present effort, the linearized aerodynamic and numerical models have been implemented 
into a corresponding first-harmonic, unsteady, flow code, called LINFLUX. At present the 
NPHASE and LINFLUX codes apply to two-dimensional flows, but extensions to three- 
dimensional flows are planned as future work. 

The three-dimensional aerodynamic and numerical formulations and the two-dimensional 
LINFLUX analysis are described in this report. Numerical results for two-dimensional un- 
steady cascade flows, excited by prescribed blade motions and prescribed aerodynamic dis- 
turbances at inlet and exit, are also provided to illustrate the capabilities of the LINFLUX 
analysis. The present capabilities for predicting unsteady subsonic flows produced by blade 
vibrations and acoustic excitations are quite good, but additional work will be required so 
that accurate predictions of unsteady transonic flows and unsteady flows excited by entropic 
and vortical gusts can be determined. Additional work will also be required to improve the 
computational efficiency of the LINFLUX analysis. 
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1. Introduction 
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Linearized potential methods are based on the assumption of an isentropic and irro- 
tational mean or steady background flow. Although these assumptions lead to significant 
computational efficiencies, they also limit the range of application of the resulting Analyses 
In particular a more general linearization is required to predict three-dimensional unsteady 
flows in which the effects of mean swirl are important, and two- and three-dimensional flows 
in which strong shocks occur. For such flows, the nonlinear Euler equations are required to 
model the nomsentropic and rotational mean flow, and linearized versions of these equations 
are required to model the unsteady perturbation. 

. . Recently, much attention had been given to the development of two-[HC93a HC93b 
KK93] and three-dimensional [HL92, HCL93] linearized Euler analyses. In these analyses, as 
in the earlier linearizations relative to potential mean flows, the linear unsteady equations are 
developed m the frequency domain and analytic far-field solutions are matched to numerical 
near-field solutions to limit the computational domain to a single blade-passage region of 

nite extent in the axial direction. Unlike the classical and potential-based linearizations, 
the linearized Euler equations are solved on a grid that deforms with the blade motion so 
that troublesome extrapolation terms in the blade surface conditions can be replaced by 
more tractable source terms in the field equations [HC93a]. Also, shock and wake effects 
are captured” within a conservative finite-volume discretization, rather than ‘'fitted” by 
imposing shock and wake jump conditions. Finally, because of the large number of unknowns 
in the discretized equation set, the first-harmonic linear equations are solved using pseudo- 
time as an iteration technique rather than by direct matrix inversion. 

Under the present effort work has been initiated to develop a three-dimensional linearized 
unsteady aerodynamic analysis for turbomachinery blade rows that accounts for the effects 
of real blade geometry, mean blade loading, mean swirl and strong shocks. It is based on the 
Euler equations of fluid motion and the high-resolution, wave-split, implicit time marching 
numerical scheme developed by Whitfield, Janus and Simpson [WJS88]. Here, the time 
dependent Euler equations are written in strong conservation form in terms of curvilinear 
computational coordinates and solved using an implicit time-marching procedure. A three- 
point, backward difference approximation, which is second-order accurate, is used to ap- 
proximate time derivatives, and a cell-centered, upwind, finite-volume discretization is used 
to approximate spatial behavior. In the spatial discretization, fluxes at the cell interfaces 
are evaluated using a flux-splitting technique in which the eigenvalues of the flux Jacobian 
matrices control the direction of spatial differencing. The nonlinear discretized equations 
are solved at each instant of time using a Newton iteration procedure, that involves an ap- 
proximate factorization of the residual equation. Gauss-Seidel sub-iterations are applied to 
reduce the approximate-factorization errors. 

The wave or flux splitting used in [WJS88] allows flow discontinuities to be captured in 
one computational cell, if the grid is well aligned with the discontinuity. In addition, since 
the splitting scheme is based on a local eigenvalue decomposition of the Euler equations, the 
mean flow can be interrogated to determine the shock locations. A change in sign of the 
eigenvalues associated with acoustic waves indicates a change from supersonic to subsonic 
flow, and therefore the existence of a shock. This feature could facilitate shock fitting in the 
future via the imposition of shock jump conditions at mean shock locations. 

Based upon the foregoing numerical scheme, Huff, Swafford and Reddy [HSR91] have con- 
structed an implicit, multi-block, finite- volume analysis and computer code, called NPHASE, 
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for time- accurate resolutions of nonlinear, U— 
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2 . Unsteady Flow through a Rotating Blade Row 

2.1 Physical Problem 

We consider time-dependent adiabatic flow, with negligible body forces, of an inviscid. 
non-heat conducting, perfect gas through a single rotating blade row such as the one shown in 
igure 1. In particular, we will be concerned with a blade row operating within a cylindrical 
duct of radius r = r D (£) and mounted on cylindrical hub of radius r = r H (£), where £ 
measures axial distance and r measures distance radially outward from the £-axis. The 
blade row consists of A distinct blades which rotate about this axis at angular velocity fi. 
In the absence of unsteady fluid dynamic forces, the blades are assumed to be identical in 
shape, identical in orientation and equally spaced around the rotor. The mean or steady 

state position of the nth blade, in the rotor frame of reference, is described by an equation 
of the form 


27T7T 

F(r, 0,£) = 0+ + /(r,£) = 0 , n = 0, 1, 2, . . . , N - 1 , (2.1) 

where 9 measures angular distance in the direction opposite to the direction of the rotation 
and n is a blade number index which decreases in the direction of rotation. 

For aeroelastic and aeroacoustic applications we are interested in a restricted class of 
unsteady flows; i.e., those in which the unsteady fluctuations can be regarded as pertur- 
bations, not necessarily small, of a background flow, that is steady in the rotor frame of 
reference. Thus, we will consider primarily situations in which the blades rotate at con- 
stant angular velocity and in which the background flows far upstream (say £ < £_) and 
far downstream (£ > £+) from the rotor are at most a small perturbation from a steady, 
axisymmetric, swirling flow. This perturbation stems from the interaction between the fluid 
and the blading and is steady in the rotor frame of reference. 

The time-dependent or unsteady fluctuations in the flow, as seen by an observer attached 
to the blade row, arise from one or more of the following sources: prescribed vibratory blade 
motions, and prescribed upstream total temperature, total pressure and static pressure dis- 
turbances and downstream static pressure disturbances, that carry energy toward the blade 
row\ For the most part, we will restrict our consideration to temporally and circumferentially 
periodic, unsteady excitations that are of small-amplitude. 

For example, we will consider prescribed blade motions of the form 

72.(r,0 + 27rn/A r ,£,f) = Re{R(r,0,£)exp[?(u;f + ncr)]}, xonfl, n = 0,l,2,... ,A'-1, 

(2.2) 

where the vector 7L is the displacement of a point on a moving blade surface relative to its 
mean position in the rotating frame; R, (|R| ~ O(e)), is a complex displacement-amplitude 
vector; us is the temporal frequency of the blade motion relative to an observer in the rotor 
frame; o is the phase angle between the motions of adjacent blades; /?e{ } denotes the real 
part of { }, and B denotes the zeroth (n = 0) or reference blade surface. 

The unsteady flow's in the far upstream and far downstream regions are in part, pre- 
scribed as a fluid dynamic excitation and, in part, depend upon the interaction between 
the fluid and the blading. Typically, an unsteady aerodynamic excitation is represented 
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as a linear combination of fundamental disturbances that are harmonic in time and in the 
circumferential direction. For example, a fundamental pressure excitation would be of the 

f ° rm p/ lT oo(x,t) = /?e{p/, T co(r)exp[i(u;t + m0) + Xm^]}, £ > & t 2 ' 3 ) 

it „ fx t) is an incident pressure disturbance, i.e., a pressure disturbance that travels 

J 1 blade row from for upstream (* < i-) or fax downstream tf > *+)’ " 15 the 

towards the blade row irom i H ' . g the number of comp l e te disturbance cycles 

temporal frequency of the dis ’ . ^ 0 < 0 < 2 7 r. The quantities w and m are 

thal Zi r r w and il ai determined from the equations that describe the fluid 

t tt'Sd Note that the interblade phase angle, o, of an incident disturbance 

i$ ^(Lelastic and aeroacoustic investigations, the applications of primary interest here 

various prescriptions of the foregoing excitations. 

Two-Dimensional Appvoximation 

, the ^ 

the physical problem just descnb . development of a three-dimensional 

Z'^r respectively; T t 

Srlt ^L^tt^ute velocity of the blade row in the 
circumferential direction, and is constant, if O 0. _ e y y etc and regard 

If we unwrap the **££££ Id tangen'tiS dLtions', respect, vely' we 

are'eff^tWely Considering the flow througharectil^ 

positions of the blade -dror ^ „ is the blade numb er index, 0 ,s the cascade 

Llg^ld G = Ge,’is the cascade gap vector which is directed along the ,-axis with 

magnitude, G, equal to the blade spacing. 

The vibratory blade motions, i.e., 


7 £(x + nG,t) = Re{R(x)exp[i(wt + nor)]}, 


G B ,n = 0,±1,±2. . . 


(2.4) 


, ... „ v aT1 j t-iinp t and for a two-dimensional flow, small 

are prescribed functions of posi i ’ ’ ted as the sum of independent 

unsteady disturbances far from the blade row can be ‘ ™ S ” , () ^ ortical , 

entropic vortical * 1 acoustic mo es^of and «, that 
sa'tTsfy the fleSTequaLs and describe disturbances that carry energy toward 
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tr;; o i q fund “ “ . can be 


and 


5-oo (x, t) - expft^ic.oo • x + w f)]} , 

C-oo(x, t) = RefC.^ exp[i(K_ TO • x + u;*)]} , 


(2-5) 

(2.6) 


(2.7) 


P/ tT oo(x, t) - Re{p/ iTco expf^Kzpoo • x + art)]} . £ $ ^ [2 yj 

ot in t dent entropic - ™ rti “‘ “ d 
'ZTJ* “ pV ft dis t 

« 2ZV7TZZ r n Zr er v an in f ententropic « 

but a more complicated relationship Exists between " md VeI ° City ’ 

disturbance [Ver89bj. 311(1 K:p00 *° r an inci< ^ ent pressure 

2.2 Mathematical Conventions 

xctssrs ~y.~z:ss7. 

(V*) 2 fk P \ ) i P cine internal energy and enthalpy with respect to 

:r :r^:: d " t r:rsti r tions - Ako ' in the 

*v >r 

fluid property, and, with the exception of the density UDDer and 1 dependent 

be used to represent the various nonlin2r^ ^ ^ 

terparts, respectively. Thus, for example, P and P wUlrefeTto t £ IZt 1 T* °T' 
fluid pressure and the pressure in a background flow that is steadv in the^ f 1 ? 6 ' ependen ^ 

the flow variables far upstream and far downstream from the blade row ° ™ ““ ° f 

- -ipt letters 

JS&TT.T “ <s T: T ‘^trsa.-srzs; 


+ tefirCeC 1 f U f fits 
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space-fixed or inertial reference franre. ^ 

rotating cylindn«ljpat.J coo. ^ma ^ measure distance along and normal, 

r, 9 A , £, where 9 - «J ' . , f , _ nffU i ar coordinates, 0 and measure angular 

respectively, to the axis of ro a ion rect i 0 n of rotation. We will also use rotating 

distances in the direction opp^rte “ h ^^ rdinates in which the or x^-axis 
(xi,x 2 ,x 3 ) and space-fixed ( 1 , 21 3 ) Q The Cartes ian (xi,x 2 ,x 3 ) and 

coincides with the axis of rotation, i.e V - lie< *1 

cylindrical (r,0,£) coordinates of a point are related by 

( 2 . 8 ) 


Xl = £ } x 2 = rsin0 and xz — —r cos 9 


or 


€ _ Xl , 9 = tan- l (-X 2 /x 3 ) and r = [(x 2 ) 2 + (x 3 ) 2 ] ( 2 - 9 ) 

We will require equations that describe the flmd V • F and 

reference frame and the rotating reference functions, respectively, are invariant 

V x F, where / and F are arbitrary scalar and „ is the total time 

under a transformation from ° n€ ' ^^Tderivatives of scalar Lnd vector quantities 

- invariant under such a transformation, 

but are related by 


dj_ 

dt 


lx A 


d l 

dt 


- (f? x r- V)/ , 


9F 

dt 


and 


d f 

dt 

d F 

dt 


— ( C2 x r ■ V)F + f? x F, 


( 2 . 10 ) 


= - + flx F 
dt 


Ut- 

, , ~,m and dFIdt are the time derivatives as seen by observers in the 

Here, for example, dV/dt\ A and a / ofrnTT .„ resoectively If the position vectors, 

and rotating frames are related by 

( 2 - 11 ) 


V A = 


Dx a 

Dt 


— + f2x^x = V + fixr, 


DV a 

Dt 


£ + nx)(V + fl*r)a + !flxV + ftxr + nx(l)xr), (2.12) 


and 


= Vx V ;1 = Vx(V + fixr)=VxV + 2fi-C + 2fi 


(2.13) 
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Here V- 4 , A a and C are the velocity, acceleration, and vorticity, respectively, as seen by an 
observer in the inertial frame; V, A and £ are the velocity, acceleration, and vorticity seen by 
an observer fixed in the rotating frame; is the angular acceleration of this rotating frame; 
and r is the radial position vector. Quantities measured with respect to the inertial frame will 
be referred to as absolute flow quantities; those measured in the rotating frame, as relative 
flow quantities. The terms 21? x V, Q x (1? x r), and 1? x r in equation (2.12) represent 
Coriolis, centripetal and Euler accelerations, respectively. Note that, if V x V 4 = 0, the 
flow is irrotational in the space-fixed system, but it is rotational with £ = — 21? relative to 
an observer in the rotating frame. 


Spatial Independent Variables 


To describe flows in which the fluid domain varies with time; e.g., unsteady flows over 
vibrating solid bodies, it is useful to consider two sets of independent variables, say (x. t) 
and (x, t), where the position vectors x and x refer to locations within the rotating frame 
of reference. The vector x(x, t) refers to the instantaneous location of a moving field point, 
say V , at time f, whereas x refers to the reference position of V. Thus, we can write 


x = x(x, t) = x + 7£(x, t ) , 


(2.14) 


where 7?.(x, t ) is the displacement of the point V from its reference position. 

This representation is motivated by the need to conveniently account for solid body 
motions in unsteady flows. In particular, the vectors xs{x s, t ) and Xs denote position vectors 
to a point on a moving surface S and to the location of this point on the reference surface 
S, respectively, and IZs = Xs — x s is the displacement of a point on the moving surface S 
relative to its stationary or reference position on S. It is assumed that TZ(x, t ), is continuous, 
single valued, and that (2.14) can be inverted to give the reference coordinates x, of a moving 
field point in terms of its current position, x, and the time f; i.e., the functions x = x(x,t) 
are also continuous and single- valued. 

The behavior of fluid properties may be followed along the path of the moving point, 
V, or at fixed spatial locations. For example, the fluid density, p, may be expressed as a 
function of x and t, where x is the reference location of V, or as a function of x and t , where 
x is the position of V at time t. Associated with these two descriptions are two derivatives 
with respect to time; d( )jdt | x , the derivative with respect to time keeping x constant, and 
d( )jdt | x , the derivative with respect to time keeping x constant. The former is the rate 
of change in a quantity as observed when following the motion of a moving point; the latter, 
the rate of change as observed at a fixed point x in space. These derivatives are related by 


d( ) 


dt 


d( ) 


dt 


+ (£-V x )( ) 


(2.15) 
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3. The Nonlinear, Time-Dependent, Aerodynamic Equations 


The description of the foregoing physical assumptions and mathematical conventions 
brings us to a convenient starting point for presenting the equations that govern the u 
fl S thrnutrh the blade row These are determined from the conservation laws for mass, 
momentum^nd^nergy. and the thermodynam.c re.at.ons for a perfect gas. We will wnte 
The governing equations, first, in terms of dependent and independent vanables measured 
with respect to an inertial reference frame, and then, in terms of vanables measured wrth 
respect to a reference frame rotating with the blade row at angular velocity fl. 

3.1 Conservation Laws and Thermodynamic Relations 

We consider an arbitrary moving control volume, V(t), which is bounded by a moving 
control surface S(xM) = 0. The conservation laws for the fluid within V at time 

written as , i . . , 

^ [ pdV\ + / p(V A - K ) ■ ndS = 0 , 

1 Jv L Js 


(3-1) 


dt 


and 


— [ p\ A dV | + [ ~p[V A -il A )}-n d S - - j s 
dt Jv [A Js 

j t pE$dv\ A + J s pEi(V A - 11*) ■ ndS = - j s 


Pn ■ dS 


P\ A • ndS 


(3-2) 


(3-3) 


dt Jv 1/1 J * 

Here dldt\ is the time derivative taken relative to an observer fixed in the inertial reference 
hame inb V J , P and E A = E + (V A f/2 are the fluid density, absolute velocity, pressure 
Id absolute specific totai internal energy, respectively; ** is the absolute velocity of a field 
, . t . tro i surface S' V A - jl is the fluid velocity relative to the moving 

££ SlkW- normal vector to the control surface; and ® denotes 

‘^Ti^firsT^Xsecrnid^ermsViTthe^'efLhand-sides of the integral mass (3.1) momentum 
„ ,T Id Irgy (3 3) conservation laws represent the time rate of increa* o a conserved 
\Z H e mass momentum or total internal energy) within the control volume and the 

tid^thmt 

the contro 'fume ’ viscosity, discontinuities in the flow variables can occur. The 

re ;^tm^"-ion laws are therefore required to describe the flow over the 
m , . r ■ n f Pres t These forms provide corresponding differential equations in regions 

t“" Jr Tondnuolly differential, e and W conditions at surfaces 
Irol which The flow variables are discontinuous, i.e„ at vortex-sheets, and at shocks. 

Rotating Frame Vanables 
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relations: 

V - 4 = v + Q x r 

E} = E + (V A fn = £r + V.J7xr+(fix r ) 2 /2 (3 ' 4) 

into equations (3.1) - (3.3). Here V is the relative fluid velocity and E T = E + V 2 /2 is the 
relative total specific internal energy. After carrying out the algebra we find that 

SiX' wv + (3.5) 


d 

It X mv + f s W ® (V - R) + PX\ ■ ndS —Ox/ (iVdV- + J? x j /ftfl 

-J s H(n x r) ® (V - R)| • ndS, 


x r)dV 


and 


(3.6) 


1 X^V + / [Mr(V - R) + PV) . ndS = -± jjjfr -Ox r + (fi x r) J /2]rfV 
-|()[V-l)xr + (Ox r) 2 /2j(V - R) + P[2 x r) ndS, 

where d/dt is the time derivative measured by an observer in the rotating frame and R 'is 

e relative velocity of a point on the moving control surface <S(x,t). Note that the left- 

hand^des of (3.5 - (3 7) are identical in form to (3.1) - (3.3). The source terms on the 

right-hand-sides of (3.6) - (3.7) account for the rotation of the blade-row fixed reference 
Irame. 

The integral conservation laws can be written in various alternative forms. For example 
we can interchange the order of time differentiation and volume integration for the volume 
integrals appearing in (3.5), (3.6) and (3.7) by employing a modified version of the Revnolds 
transport theorem (see [Ari62] p. 85 or ([Mou84] p. 51), i. e ., 


d_ 

dt 


X< )dV -L 


'd( ) 


dt 


+ v x • [iz ® ( )] dv , 


(3.8) 


vvhere h <g> ( ) reduces to 7L( ) if the quantity within the parantheses is a scalar. We can 

also apply Green’s theorem, i.e., 


j s n • ( )dS = V x • ( )dV 
to convert surface integrals in (3.5)-(3.7) to volume integrals. 


(3.9) 
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Thermodynamic Relations 

In addition to the conservation equations, dynlmt problem and to 

of uSady flow processes. In particular, we require the equation 
of state for a thermally perfect gas, 

and the relation between the specific internal energy and the temperature for a calorically 

perfect gas, £ = f = ( 7 _ 1 )~'p/p . (3 ‘ U) 

. anrl „ ; s t he specific heat ratio (constant pressure to constant 
Here f is the tempera ure ,qn ^ (3 that the pressure can be expressed 

volume) of the fluid It follows from (3.10) and (3jl) P 

in terms of the dependent variables />, V and b T />, _ 

p = ( 7 - V)p\Et - (V'' 4 ) 2 /2] = (7 - 1 )KEt - V / 2 ) l ' ’ ~ 

variables p, V' J » and Of' or P, V, and _£r- 

The enthalpy, A, of a fluid particle is defined by 


H = E + P/p , 

and it follows from (3.10), and (3.11) 

H = T = 7 _1 (7 - 


(3.13) 


(3.14) 


Although enthalpy and temperature are ifferent^md ^“™ me ^7 nali2e , he valio us flow 

i: r » - - « - - * rn5) 

H* = H + (V A ?I2 = % + PIP 11 

and /- £rA v . o v r - (Q x r ) 2 / 2 . (3»16) 


fl x =ff + CV2 = & + P/P=HT-' , -nxr-(nxr) 2 /2. 

We will also find it useful to introduce the fundamental thermodynamic identities 

--1 j n 


f dS = dE + Pd{p 1 ) = dH — p dP, 


(3.17) 


„ 4 f fll ,id It t hen follows from (3.10) through (3.17) that 

tTnlpyt^rth^sure and density by the differential expression 


dS = 7 -1 dP/P - dp/p • 


(3.18) 


P t - n IS) can be integrated to relate the change in entropy between ay 
Equation (3.18) can b 6 changes in pressure and density; i.e., 
rium states to the corresponding changes in p 

Pp~^e~^ § ~ {P'p ' 1 e" iS )nei , 

where the subscript refers to a reference thermodynamic state. 
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3.2 Field Equations 

Equa tio n s that describe inviscid fluid motion at the field points x - S + WY* « • 
“'rr 8 ' “« fro™ ‘he integral conservation lawL e 7 (3 1 W 3 ' f 

the °t r of ^ 

integrals. After taking the limit of the voluie IntepS asVllf^ '" tegralsto volume 
conservation laws, expressed in terms n f th a F V ^' we find that the 

independent variables! and tk^Zlo^ “ * V ^ & and the 


dp 

dt 


dp V 

dt 

and 

dQEr) 


+ v x • (**) + V x • [(V - n)p] = 0 . (3.20) 

+V > .(W® J 5V) + V x .[(V-R)®p-V + p II = -^OxV+Oxr+OxfOxr)] (3.21) 


dt 


+ V, 


■(npEr) + V x .((V--k)-pE T +PV] =--p[v.(nxr)-n 


x r • i? X V 


The right-hand-sides of the momentum (3.21) and enercv f3 221 *• , , ( 3 - 22 ) 

phfied by making use of (3 20) and f3 201 ^ ‘ 22 equatlons have been sim- 

algebra. The source term on the “IS H ^ performing some 

for the Coriolis, Euler and centipefdtcTele^ amentum equation (3.21) accounts 
vide a complete set of govern! ^tionf X ^ partide ' T ° P"> 

the thermodynamic relation (UVwhT^ZTT ” ^ ^ by 

total specific internal energy id ““ *"** "■*« 

TtXicemlnt ££ #£“**“* «“ fl “ id 

a 

(3.20) - (3.22). However, at^s hat cZn to ita”^ T ^ from 

a direct correspondence between the differential andThe inteii [(H) « 711 ?" “ “TJ* 

conservation equations. Also, equations (3 201 - f3 221 h a k 1(3-5 • (3 ' ? ^ f ° rmS ° f the 
independent variables x and f Til ,11 1 7 j , (3 ‘ 22) have been wntten in terms of the 
independent variables by introducing theroatM aWS could ^ be written using x and t as 
adopt this strategy in ^ *' ^ 

merical simulations so th7t Mdilfi^teTtpTured”' but toe Preferred l °' 

ThtTouTfcr^o^^kft^p'r 15 ' in m * Mkua * ^ “‘ k toSstssi! 

fve form of the mis balaie equli T ^ "" St ™* b ‘“ d ^a. The convec- 


£l 

Dt 


7T- + ^V x • V = o 


(3.23) 
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(3.24) 


, _ _ g/m + V • V, is the material or convective derivative operator relative to 

" vfr fTxe/fn tht rotating frame. The convective form of the momentum equation, 

+ V X P = -p [212 x V + 1? x r+ « x (« x r)] , 

,s obtained by combining (3.21) and (3.23), and the convective form of the energy equat.on, 

ie " iM = -.PV„.V, < 3 - 25 > 

' Dt , 

is obtained by combining (3.22), (3.23) and (3^24) jmd the 

Alternative convective forms^ the energy q ^ themodynamic relations (3.17). 

thermal energy equation (3.25) with V ) 

We find that - _ D § jjft DP _ n (3.26) 

^ + PV x V = pT- 6r = P- 5r -- 5r -°- 

r u a„\A nartirle must remain constant in continuous regions 
In particular, the entropy of each fl P , , t ■ h regions, the pressure 

ihsss s £ 's:^Uon -u >r - 

constant. 

Surface and Far-Field Conditions 

For application to turbomachinery blad^rows 

plemented by boundary conditions a - f t ady aerodynamic behavior is usually 

at the inflow and outflow boundaries. Trans ent unst«dy Therefore, for 

not of interest in turbomachinery and i$ „ ot requ ired. 

SUC F„Tfnr“d n flowsTefowUngency condition appUes at the moving blade surfaces, i.e., 

(3.27) 


(Y _ fz,) • n = 0 , x £ B n ( or * ^ ’ 


and at the stationary duct walls, i.e., 


V . n = 0 , at r = rtf and r — rp- 


(3.28) 


The blade velocity H in (3.27) is -P^ribed quantity l<* (£) from the blade 
We also require conditions on e of the compu tational domain. Typically the 

row i.e., at the inflow and outflow bound f th- total pressure total temperature 

circumferentially- and temporally-avera^e ^ ^ s ° 0 f ra dius at the inflow boundary. At 
and the inlet flow angle are specified ,, a * d temp orally-averaged pressure is specified, 
the outflow boundary, the circum er ^ pressure and total temperature flue- 

consistent with radial equilibrium. ^ and exit , that carry energy towards 

tuations at inlet, and the pressure should be noted that the unsteady excitations 

the biade row, must be specified As discussed in §5, such 

— - - - 
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and far downstream of the blade row. Total pressure and total temperature fluctuations at 
exit and unsteady pressure disturbances at inlet and exit that carry energy awav from the 
blade row must be determined as part of the unsteady solution. This is usually accomplished 
by extrapolating numerical “near-field" solutions through the inflow and outflow boundaries 

shod.-s P K n . C Th e ’ JUm ? COnd, * lons shcmld be 1 txi posed at vortex sheets (“viscous lavers") and 
cks, but the usual procedure in numerical calculations is to capture such discontinuities 
by solving the conservation forms of the field equations (3.20) - (3.22) over the entire fluid 
domain. However for the sake of completeness, we will describe the jump co"d,W that I 

laws 3 1) (3 3)' W3 -n h b Se C0 " d d tiOnS ^ 0bUined fr ° m the mtegral 

(51) (3.3) or (3.0) - (3.1) by considering a control volume that contains a segment 

such a surface, and taking limits, first, as the lateral extent of this volume normal to the 
surface segment approaches zero, and, then, as the area of the segment approaches zero 

acr JsTuihXT u S f ° r T 18 ™* mass - and energy at a surface. S. 

across which the flow variables are discontinuous are 


and 


[tf/I = o , 

M f [V A ] + [PJn = M/fVJ + [Pin = 0 


(3.29) 

(3.30) 


+ [PV^J • n = M s {E t } + [PVJ . n = 0 , (3 31) 

respectively. Here [ ] denotes the jump in a flow quantity as experienced bv an observer 
when moving across the surface of discontinuity 5 in the n-direction and 

Mj = ~p(V A - TZ S ) ■ n = p(V - il s ) • n , xeS (3.32) 

is the fluid mass flux through the surface. 

Equations (3.29) - (3.31) have been derived by setting the field point velocity 7Z in the 
in egral conservation laws equal to the surface velocity 7l s . Thus, these equations are applied 

I t P .° in S ’ X f ^ X . 6 S ^’ that Iie on a movin g “boundary layer”, wake or shock surface 

at ^TT p PP1Cati r the SUrfaCe VdOCity Vect ° r is P rescribed at blade surfaces, but 
at wake and shock surfaces it must be determined as part of the overall time-dependent 
solution for the unsteady flow. p ’ 

the condit™ r ns eX ^ SUPP ° r ‘ * * tan « ential velod ‘y. !•«•. WJ • r # 0, it follows that 

Mj = XV - 7l s ) • n = 0 (3.33) 

I-Pl = 0 (3. 34 ) 

and 

[VJ • n = 0 , (3.35) 

prevail at vortex sheet boundary-layer and wake surfaces. Since the vortex-sheet boundarv 
lavers coincide with the blade surfaces, (3.33) is equivalent to the inviscid flow tangency 

from H 30Wb f ^ ^ 1$ generaIly nonzero (be-, M, # 0). Hence, it follows 

across ' ho S com P one ^t of the fluid velocity V • r, must be continuous 

across shocks. The remaining jump conditions, along with the thermodvnamic equations of 

state, are hen required to determine the shock velocity, and the changes in the normal 
hr3rt4°shoc e k 7 and the therm ° dynamiC Pr ° PertieS ° f the fluid as ^ Passes 
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3.3 Matrix Forms of the Governing Equations 

The conservation equations are ^often “^““J^Ling differential expressions 

the integral conservation laws (3.5) - (3.<) 

(3 ' 20) (3 '" ) a £ f VdV+ (W,- VKtJn^dS = j v SdV 

dt Jv J s 


(3.36) 


and 


dV 

dt 




(3.37; 


tt t, a „j C are referred to as state, extended flux, 
respectively. Here, the column vect ° rs ' j’_ 9 3 are the Cartesian components of the 

and source term vectors, respectively, K Xj ,J , 


vector 1Z, and a summation over r ^P e ^ ted mdl ^ S p ■ ^ \ 2 3 and the source term vector, 
The sUte vector, U, the extended flux vector, F„j - 

S are defined by 


and 



'p 

• 

l 

*-> 


pV xl 


}V„V,, + 

u = 

pV X2 

, ij = 

pVx 3 Vx, + 


pV X3 


pV„Vx, + PS 3, 


L P&T 

. 

[pEx + P)Vxj - 

- 


0 




0 


s = 


p{2nv X3 + fa* + ft 2 ^) 



— p(2flVx 2 + tlx 2 - H 2 x 3 ) 


’ _ 1/2 ... , 5, and the source term 


v I2 - x 2 Vz 3 ) + « 2 (^v ; 2 + *&>)] 

cpressed as explicit functii 
vector S can be expressed 


(3.38) 


p[n(x 3 i 

The flux vectors F, can also he expressed as * ^“ttn of the 
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state variables and the spatial coordinates x,, i = 1,2,3, i.e 



0 

0 

0 . 2 (j]X 2 + 2 O.C 4 + flUxXs 

n 2 6ix 3 — 2 fib' 3 — nu lX2 

n 2 (u 3 x 2 + u 4 x 3 ) + fi(6' 3 x 3 - u 4 x 2 ) 

(3.39) 


P = h- 1)[£/ S - ur\u 2 + ill + ill)/ 2] (3.40) 

Our derivation of the field equations and surface conditions that govern the nonlinear 
time-dependent, inviscid flow through a blade row operating within a cylindrical duct is 
now complete. The nonlinear, unsteady flow problem is formidable, and, to date, there 
as een relatively little attempt to apply time-accurate numerical simulations of unsteady 
flows m turbomachinery aeroelastic and aeroacoustic response studies, not only because of 
the assodated mathematica 1 complexities, but primarily because the computing resources 
needed for such simulations prohibits their use in detailed parametric investigations. Thus, 
in the next section we will make use of additional simplifying assumptions with the inten- 
tion of providing a useful unsteady fluid dynamic model for turbomachinery aeroelastic and 
aeioacoustic design applications. In particular, we will restrict our consideration to small- 
amphtude, temporally and circumferentially periodic, unsteady excitations, and unsteady 
perturbations of a nonlinear background flow that is steady in the blade-row fixed, reference 
rame In particular, the unsteady perturbations in the flow are assumed to be driven bv 

prescribed blade motions of the form (2.2) and/or external aerodynamic excitations of the 
torm (2.3). 
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u The Small Unsteady-Disturbance Approximation 


4.1 Mathematical Preliminaries 


tntinnal resources required to simulate nonlinear inviscid and viscid unsteady 
The compu rrm tmue to prohibit the use of such simulations in repetitive 

fluid dynamic behav.or will continue P approximate analyses, based on the 

zzs+ssz « r“u d^:“— 

““ - * - — - - 

aeroelastic and aeroacoustic design calculations^ three-dimensional blade row rotating 

Thus. we ^noniinear, time-dependent, 

at constant angular velocity, \L V r- f Qt t \ << i] unsteady exci- 

governing equations, which « P" ° £ th ' *** “ 

tations. If the excitations are of p • non li ne ar zeroth-order, background 

- r “i^nce ,, 

with the blade row. • that describe the inviscid unsteady perturbation we 

To determine approximate equation acv mntotic series in e Thus, we consider 

must first expand the dependent flow variables m asymptotic series 
dependent variable expansions of the form 

V[x(x,i),t) = V(x) + qx(x,t).t] + ...=V(x) + vlx(x,«Ml+... . (4-1) 

Here V(x,() is the nonknear "“ te ^ t ^“^or skLdTvelocity at Hie reference position 

X ofC pi? y[x(x,i),t) ~ OM Vf terms, 

velocity at x at the mstantanMus^ ^ or steady . sUte position, and assume 

IS “It A av/atu = 0 and therefore, it foiiows Horn 

(-2,15) that ^ \ becomes identical to the asymptotic expansion used 

If we set n = 0, then x - x, i > ^ linearizations relative to nonuniform 

^^Zs “eXvt 9 3, : The 


p J = (P + 'R-VF + ...) Is 


(4.2) 


are reqmred to refer iuformat ion at moving «ade surfaces to the Unde ^ition. 
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surface. These terms are difficult to evaluate accurately. Consequently, recent linearized 
sol ut ,o„ methods [HC93a, HC93b, KK93. HL92] have employed thelore gener If 
expansion (4.1). With this expansion, if the displacement field TC(x,() satifies the condition 
/ _ K Bn {x t) for x € B nt extrapolation terms will not appear in the blade boundary 
conditions or surface pressure expressions. Instead, source terms, which are easier to evaluate 
accurately, will appear in the field equations for the linearized unsteady flow. 

t'p to this point we have expressed the fluid dynamic equations in terms of the indepen- 
dent variables x and t. However, the use of dependent variable expansions of the form (4 1 ) 
renders it more appropriate to express the steady and linearized unsteady equations in terms 
of the independent variables x and t. The temporal derivatives d( )/dt\ and d( )/dt\ 
transform according to (2.15); the gradient operators, according to 


V x ( )«V*( )-(V,< 8 >rc).V x ( ) + .... (4.3) 

The equations that govern the zeroth-order steady and the first-order unsteady flows are 
obtained by substituting the series expansions for the dependent flow variables "e g (4 1 ) 
and the independent-variable relations (2.15) and (4.3) into the nonlinear, time-dependent, 
governing equations; equating terms of like power in e; and neglecting terms of second and 
higher order in e. After carrying out this procedure, we find that nonlinear and linear 
variable-coefficient equations are obtained, respectively, for the zeroth- and first-order flows. 

e variable coefficients that appear in the linearized unsteady problem depend upon the 
zeroth-order, steady background flow. 

We will write the steady and linearized unsteady equations in terms of the rotating frame 
dependent variables, i.e., the relative velocity [cf. (4.1)] and the relative total specific internal 


-SWx, t) - E + V 2 /2 - E + V 2 {2 + e + V • v + . . . = E T (x) + e T (x, t) + . . . . ( 4 . 4 ) 

It follows from (3.4) that the perturbation absolute velocity and specific total internal energy 
are given by 6 " 

v- 4 (x, t) = v(x, t) + O x (r - r) 

(4.5) 


24 - 


= e + \ A •v i, = ej' + v- flxr + V- ftx(r-f), 


where r = r(x) and r = r(x). 


Consequences of the Small Unsteady Disturbance Approximation 

A significant simplification offered by an unsteady aerodynamic linearization is that the 
fluid dynamic responses arising from various sources of unsteady excitation are not coupled, 
and hence, they can be determined separately. Indeed, it is usually sufficient in turboma- 
chinery applications to develop solution procedures for a single harmonic (in space and time) 
component of a given unsteady excitation. Solutions for an arbitrary independent excitation 
and lor arbitrary combinations of various independent excitations can then be obtained bv 
superposition. It should be noted, however, that the linearization does limit the unsteady 
aerodynamic response phenomena that can be analyzed, since nonlinear inviscid unsteady 
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r„« d ofToc^ cannot be taken mt0 

account. . . , the rotor fixed frame and the equations that govern 

Since the background flow is steady fluctuations in the unsteady flow prop- 

the first-order unsteady flow in tim e must have harmonic time 

* — . a — «— ■ 

for the first-orcler flow properties, e.g., ^ 

v(x,f) = He[v(x)ex p(wt)], 

where * is the temporal 

and adopting the convention a The complex representation (4.6) can e us f 

«- "*■ 
facilitating the determination of * rotot geometry and the mom mlet and 

for example, we can write _ _ _ / 4 7 ) 

V(f, 9 + 2ir nfN, 0 = v ( r ’ ^ & ’ 

where J. f are the cylindrical coordinates t 

— — — -*•’ (4 8) 

V(f , 0 + 27rn/dV, () = v(f , 9 , 0 exp(zna) . 

Here . = 2 *m/N is the phase angle « 

sages, N is the number of distinct b a es m Conditions (4.7) and (4.8) allow a 

of disturbance cycles over one revolution 0^_ flows tQ be limited to a single 

numerical resolution of the steady a properties at a point in the nth blade passage 

- - * — 

(n = 0) passage. 


4.2 The Steady Background Flow 


4.2 lne oieau y — ©- - , i 

The field equations and surface con<Uti«M' that^govern time- 

ground flow are obtained by substituting ^ zeroth . or der terms. The zeroth-order 

dependent nonlinear equations and re am.ng V turbation . Therefore, a solution to 

steady background flow are _ ^y) _ 0 , ^ 4 ' 9) 
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and 


Vx • (V ® pV) + VxP — —p[ 2f2 xV + fix(Jlxr)l, 


(4.10) 


• (VpEr + PV) = p£l x r • f2 x V . (4.11) 

The convective forms of these equations are 

V • Vx/j + p'Vx ■ V = 0 , (4.12) 

p(V ■ V*)V + V*P = —p[2Cl x V + x (fi x r )] , (4-13) 

and 

V-Vx5 = 0. (4.14) 

Here p, V, P, E T and 5 are the density, relative velocity, pressure, relative total specific 
internal energy, and entropy, respectively, in the steady background flow, f = r(x), and the 
steady pressure can be determined in terms of the dependent variables p and E T or p and 
5, by using (3.12) or (3.19) respectively. 

Surface and Far-Field Conditions 

In principle, surface conditions for the steady background flow follow from (3.29) through 
(3.35) and are imposed at the mean or steady-state positions of the blade, wake and shock 
surfaces. Mean blade positions are prescribed, but the mean wake and shock locations must 
be determined as part of the steady flow solution. However, the usual practice in nonlinear 
inviscid calculations is to capture shocks and wakes by solving the differential forms of the 
conservation laws throughout the fluid domain. Since we will follow this practice in the 
present study, shock and wake conditions will not be explicitly imposed. 

By assumption, the steady background flow remains attached to the blades and the duct 
walls. Therefore, the conditions 

V • n = 0, for x 6 B n , r = r H , and r = r D (4-15) 

apply at the mean blade surfaces and at the inner and outer walls of the duct. 

Conditions must also be imposed on the steady background flow far upstream and far 
downstream from the blade row. Steady-state non-reflecting inflow and outflow boundary 
conditions can be constructed [SG91] in which circumferentially averaged values of total 
temperature, total pressure and the flow angles are specified as "functions of radius at the 
inflow boundary, and the circumferentially averaged value of the pressure is specified at the 
outflow boundary, consistent with radial equilibrium. Circumferential harmonics of the mean 
flow variables are not set to zero but are allowed to evolve to values that are consistent with 

the existence of an infinite annular duct upstream and downstream of the computational 
domain. 

Finally, the periodic condition (4.7) can be imposed to restrict the steady flow solution 
domain to a single extended blade passage region, i.e., a region of angular pitch A 6 = 2ir/N. 
Since inlet and exit conditions are imposed at finite distances from the blade row, say at 

£ ~ an< ^ £ = £+> respectively, the solution domain is a single extended blade passage 
region of finite axial extent. 
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4.3 The Linearized Unsteady Flow 

The differential equations that govern the first-order or linearized inviscid unsteady flow 

, .1 hv substituting the asvmptotic expansions for the flow variables (e.g., ( • )) 

are determined by substitu g - the full time-dependent governing 

and the independent variable transformation (x ) 90W3 91) 

steady disturbance field: 

(4.16) 


dp 

dt 


+ Vx • (pV + pv) = ii • VxP + [(V* ® ^) ■ • (p v ) 0 


and 


—{pV + pv) + Vx • [V ® (pV + pv) + V ® pV + pi] 

lx 

= n • Vx(pV) + [(Vx ® K) • Vx] • (pV ® V + PI) 
-p[2U xV + flx(flxr)]- p[2fi x v + Cl x (fi x Ar)] 

—{pEr + per)! + Vx • [V(p£x + per) + *p£r + P * + P Y 1 

lx 

= n . Vx(pPr) + [(Vx ®^)'Vx]- [{pEr + W 


(4.17) 


(4.18) 


+p[flxf-flxv + flxAr-flxV] + pflxr-flxV. 

Here p. *. P and .r are the time-dependent 

and relative total specific thlf’the first-order pressure, density and 
expanding the equation of state (3.14), we 

total internal energy are related by 

p = ( 7 - l)[p(e T - V • v) + p(£r - ^V 2 )] • < ' 4 ' 19 ^ 

JeglrdeTitur ctTerms andtave'ttXTe’en oTthe rightW- sides of equations 

( 4 . 16 ) — ( 4 . 18 ). 

Convective Forms of the Linearized Unsteady Equations 

( £■ + nr dpr unsteady equations follow after applying the 

The convective forms of the ■. ( ^ 1 wiab L an d the independent variable trans- 

asymptotic expanses for the detent no We find ^ 

formation x-^xto equations (3.23), (3.44) ana to ) 


Dp 

Dt 


+ ».v,? + AV,.» + ^.-v-*-Vrf + [(v,®ii).v, MW, 


(4.20) 
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p 


Dv 

^- + (v - Vs)V 


+ MV - v*)v + 


and 


= p{ K ■ Vx)V + p[V • ( V* 0 K) • V*]V + (■ V* 0 7L) • V X P 
-/?[2Q x v + ft x (ft x Ar)] - p[2ft x V + ft x (ft x f)] , 

Ds 

— + (v • v,)5 = n ■ Vx5 + [V • (V, 0 n) • v*]s- 


(4.22) 


where 5 is the first-order entropy. In equations (4.20) - (4.22) D/Dt = d/dt\ x V v is 
a convective derivative operator based upon the mean-flow, relative velocity. * ' 

t , 1° com P let ^ the s y s ^m of convection equations we require an additional equation relating 
the hrst-order density, pressure and entropy. This is obtained by expanding the thermodv 
namic relation (3.19), using (4.1), to obtain 


P = p [ 7 1 P/P -s} = A~ 2 p -ps = 0. 
Harmonic Unsteady Excitation 


(4.23) 


If we assume that the excitation driving the linearized unsteady flow is harmonic in time 
with temporal frequency u>, then we can express each unsteady flow variable as the real part 
of the product of a complex amplitude, which depends only on the spatial coordinates x. and 
the exponent, ai function e xp (z^0, cf. (4.6). In this case the local, d/dt\^ and convective. 

IDt, time derivatives of, say, the first-order unsteady pressure are given by 


dp 

dt 


X 


Re{iupexp(iut)} 


(4.24) 


and 

Dp ( f) 1 

— = Re{(iu + V • V x )p exp(fu;£)} = Re < exp(lu;f) l , (4-25) 

where D^/Dt = iu + V . V*. 

If we now replace the first-order unsteady flow properties in the linearized field equations 
t eir complex representations, we arrive at time-independent, linear, partial differential 
equations for the complex amplitudes of the first harmonic unsteady flow variables For 
example, the conservative forms of the first-order equations have the form 


tup + V* • (pV + pv) = iuR ■ V*p + [( V* 0 R) • Vx] • pV , 


(4.26) 


iu(pV -f pv) + Vx • [V ® (pV + pv) + v 0 pV + pi] 

= iu/R • Vx(pV) + [(Vx 0 R) • Vx](pV 0 V + PI) 
-p[2Ul x V + ft X (ft X r )] - p[2ft x v + ft x (ft x Ar)] , 
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and 


where 


iu{pEr + per) + V* • [V(^t + per) + v P £ x] + ^x • + pV) 

= iujR • Vx(p£x) + [(Vx ® R) • V*] • (pE t + P)V 

+p[n, X r-fixv + nxArfixV] + /)flxr-flxV 

P = (7 - l)[p(cr - V • V) + p(£r - ^ 2 /2)] • 

Boundary Conditions 


(4.28) 


(4.29) 


Conditions on the unsteady perturbation at the befwfen 

are obtained by substituting reference surface 

the unit normals at points x — x + /v^x, ) 

positions, i.e., , (4.30) 

n Bn = n B „ - Vr®'R’ n Bn + •• , 

into the full time-dependent surface condos. ; ^t-t,ngout the 

" bude surface ' Thus ’ if r and " are orthosonal 
unit tangent vectors at a point on B n , we can write 


d d 

^T = r- + u- 


(4.31) 


aS v-n = [fc + V- V T ®^] -n ,x€ B n . ( 4 - 32 ) 

Condition (4.32) indicates that the normal component 

*- th - ^ are 

assumed to be stationary, the appropriate conditions are 

v • n = 0, for r = r H and r = r D i 4 ' 33 ) 

if the unsteady excitation is harmonic in time, then the linearized flow tangency condi- 
tions can be written as 

v • n = [iu>R + V • Vj- <8> R] • n , x on B m , ( 4 - 34 ) 

and (4.35) 

an v • n = 0, r = r H ,r D v 

where v(x) and R(x), * € B, are the complex amplitudes of the unsteady velocity fluctuation 

and the unsteady blade displacement, tespect.vdy. periodicity [cf. (4.10)1 and 

linearized unsteady flow. The latter condit.ons 

will be discussed in some detail in §5. 
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4.4 Matrix Forms of the Steady and Linearized Unsteady Equations 

The nonlinear, time-independent Euler equations (3.20)-(3.22) have been written in ma- 
trix orm in (3.37). For small-amplitude unsteady perturbations of a nonlinear mean or 

steady background flow, the state vector U in (3.37) can be approximated bv the asvmp- 
totic senes - r 

U[x(x,f)] = U(x) + u[x(x,f),t] + ... . ( 4 . 36 ) 

where U(x ) is the vector of conservation variables for the steady background flow at x 
and u[x(x *_M is the vector of conservation variables for the first-order unsteady flow at 
x - x 4- 7c(x, t ). The flux Fj and source term, S, vectors can be approximated using Tavlor 
series expansions about the mean flow, U, and the reference spatial location, x. i.e.. 


and 


F ; (U) = F j (U) + ^u + ... 

S(U.X) = S(U,x) + |iu + (7e- V*)S + ... 


(4.37) 


The Jacobian matrices dFj/d U = {dF tj /dU k } and dS/dV = { dS t /dU k }, where the sub- 
scripts i and k refer to the zth row and kth column, respectively, depend only upon the mean 
flow solution, U and the reference spatial location x. 

In terms of the physical variables 


p 


PVxj + pV X] 

pv Xl + pV x i 


pv Xl V Xj 4- pV Xl V XJ + pV Xl V Xj + pSj ! 

PVx 2 + pV X2 

dFj . 

’ <9U U_ 

PVxjV Xj + pV X2 V Xj + pV X2 V Xj + pSj2 

pv X3 + pV X3 


V Xj 4- pV X3 v Xj 4- pV X3 V Xj + pS ]3 

. per + pE T . 


phrV X] 4- pH T v Xj 4- pH T V Xj , 



0 


0 


0 


0 

ds . 

du u “ 

2Q(pv X3 +pV X3 ) + n 2 px 2 

, (7£-V*)S = 

pn 2 n X2 


-2 fl{pv X2 + pV X2 ) + n 2 px 3 


pn 2 n X3 


Q 2 p{x 2'V X2 4- x 3 v Xi ) 

. +ti 2 p{x 2 V X2 +x 3 V X3 ) 


.P&(v X 2 n X2 + v X3 n X3 )_ 


where p is given by (4.19), H T is the relative total 
background flow and 


specific enthalpy [cf. (3.16)] in the steady 


h T = e T + p l (Pp/p + p) 


(4.39) 
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is the first-order relative total specific en P transformation relations 

,, r « ** *. ^ 

background flow and the first-order unsteady perturbation are governed by the equations 


dx, Fj S 


(4.40) 


and 


du 

dt 




. dxj 

x J 

lw u ) 


(<m 


• v *) Fj + fu a + ' V * )S ' 


(4.41' 


(4.42) 


ltfspctuvciv. AJ. ~ 

temporal frequency u>, we can set 

u = Re{uexp(zu>t)} and ft = Re{Rexp(tu/f)} , 

, , nt _ nf the state vector u are the complex amplitudes of the linearized 

where the components o linearized unsteady equations 

unsteady conservation variables, in this case we can 

in the form 

0 (dFj \ . — ^ tir , ( 9 J*l . T7-A p -i- — u + (R • Vx)S , (4-43) 


- - + w . [mJ a ) “ ! " (R ' V * )U + (i ' Vt ) F] + + <R ' Vs)S ’ 


The general 

aFi/wTndaS/aU that appear in the linearized unsteady equation, (4.43). 

Forms of the Equations Used in Numerical Calculations 

For numer i°^ c ^ Cu ^^ M j)^^^^g^OTrnpl > isl^l by^rmlti^lying'the'te^Mthe'^eady 

m strong the reS ult to the right-hand-side of the linearized unsteady 

^tCn (4 43) ! « «« algebra. This procedure results in a first harmonic 

unsteady equation of the form 

f dFj x 9 


*-» + ± §“ =-MV s .R)U + .^u) 


(4.44) 


+ 


_d 

dx 


d ^F m - (Vx • R)F; 


i L 


dx 


m 


d x 3S 
+ dr : {RxiS) + 


The first term on the right-hand-s.de of (4.44) is ler^ 

dependent variables within a moving contro v ° j j t [, e l as t term, with the 

-r:^hi reference frame, caused by 

the unsteady perturbation,. 
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For the numerical integration of equations (4.40) and (4.44), it is convenient to make the 
state vectors U and u pseudo-time dependent, as first suggested by Ni and Sisto [NS 761. 
Thus, we set 

U(x, t, t) = U(x, r) + i?e{u(x, r) exp(«u;f )} , (4.4.5) 

and write the steady and linearized unsteady equations in the form 

dv dFj 

JP t + aiJ = s < 4 - 46 > 

and 

du , ■ d (d Fj \ d 

F , + ’"’ u + aij ( au ' U J = • RiU + u> 

aij, 1 a as (4 ' 4,) 

^F m -(V s .R,F,j + _(^S ) + §u 

Now, both the nonlinear mean flow equation and the linearized Euler equation contain 
explicit pseudo-time derivative terms, that allow iterative solutions to be determined using 
conventional time-marching algorithms. The equations are marched forward in pseudo time 
until the steady and the complex amplitudes of the unsteady conservation variables reach 
steady state values. Because we are only interested in the steady state values of U and u, 
there is no need to march the equations time accurately. Therefore, acceleration techniques 
such as multiple-grid and local time stepping can be applied to greatly reduce the computing 
time needed to solve an unsteady flow problem. 

For finite volume calculations integral forms of (4.46) and (4.47) are required. These may 
be found by direct integration, and are given by 

^ L vdv + Is* 1 *-* 5 = l sdv ■ «•«» 

and 


d 



a + “) x ^ + x ( “) n *’ ds = X (v * ■ r)iwv 


+ Js lUjRl]V + 8x1*™ ( V ^‘ R ) F J ri X] dS+j^ —(R Xj S) + —u dV , 


where all terms on the right hand side of (4.49), except the last, depend only on the mean 
flow and displacement field, R(x). 
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5. Far-Field Behavior 


Conditions on the linearized LVe Zl'X 

outflow ({ = « + ) boundaries are needed to altow the c ^ ^ ^ prescribe incoming 

one of finite extent in the axial irec • ’ { ttie „ overn in6 equations, and the 

aerodynamic disturbances ( e * clta u^diics must be transparent to all outgoing waves, so 
computational inflow and outflow generated An approach, that has been applied 

that spurious, nonphysical reflections are « * determine analytic solutions 

successfully in two-dimensional of a blade row, and to 

for the unsteady flow variables P fi l J so l u ti on [Ver89b]. Unfortunately, analytic 

match these solutions to a 

solutions fo ^;*”'“Xro” S nd flow. Thus, at present, approximate far-field cond.t.ons 

in computational Ui— oUh^— 1 flows ^ ^ ^ 

In this section, equations that descnb y cylindrical duct will be examined. We 

field (i.e., for fe) * a blade row but as a convenience, 
will write these equations or a sta ™ < bsolute vec tor quantities and the subscript 

we will omit the superscript A m descn g Thus, in this section, 

tzz&zX « tsr^ o 

s=T t “■ - * 3 0 the 

only on radial position, the mean flow is isentrop.c, a 
steady velocity is negligible, i.e. 


V ft! V${r)eg + V^{ r ) e i ■ 


Under these conditions [cf. (4.13) and (3-18)] 


-i dP 


— = 7 ( 7 - 1 ) 

K dr n 

As a consequence of the foregoing mean 
motion (4.22), (4.21), and (4.20) reduce to 


-l d{p l P)_ _ r -iy 2 


(5.1) 


(5.2) 


dr 


-flow assumptions, the first-order equations of 


£i = o. 

Dt 


°l+v tv+m-v,) 

Dt dr 


- + r'Vx P - (prA 2 )-'V?f*r -0. 


and 


£i + pA 2 V x ■ V + r~ X pVgVr = 0 • 


(5.3) 

(5.4) 

(5.5) 


Dt 


r , . (a no\ orir) the mean-flow lsentropic 

" ? P %Ze <M> continuity (5-5) equations, and 
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D /Dt = d /dt | x + r-'Ved /dd + V^d /dS,. (5.6) 

The linearized unsteady flow in the far-field can be determined by solving the system of 
convection equations (5.3), (5.4) and (5.5), for the first-order unsteady entropy (£), velocity 
(v) and pressure ( p ). It follows from (5.3) that an entropic disturbance is convected by the 
mean flow, and it therefore has a general solution of the form 

5(x,t) = s(r,r6 - Vgt,S - V^t) (5.7) 

If the mean flow is uniform, i.e., V = V^e^, the unsteady equations (5.3)-(5.5) have 
analytic solutions, which will be examined in detail in §5.1. These give some insight into 
three-dimensional unsteady flow behavior within a cylindrical duct. Analytical approaches 
lor three-dimensional mean flows with axial shear and swirl will be discussed in §5.2. Another 
important situation in which analytic far-field solutions can be determined, is one in which 
the radial component of the unsteady velocity can be regarded as negligible [i.e., v T ~ 0(e 2 )]. 
This "two-dimensional”, unsteady-flow approximation will be examined in §5.3. Finally, the 
techniques used for formulating and implementing numerical far-field boundary conditions 
will be discussed in §5.4. 


5.1 Uniform Mean Flow in a Cylindrical Duct 


First, we consider the unsteady perturbation of a uniform flow in a cylindrical duct, with 
velocity V = V^. In this case the unsteady equations (5.4) and (5.5) reduce to 



£+r‘vrf-°. 

(5.8) 

and 

«Am 2 V x .v = 0, 

(5.9) 

where 

D/Dt = d /dt\ x + V'd /dt, . 

(5.10) 

After taking the curl of each term in the momentum equation (5.8), we find that 




(5.11) 


where C = V x v is the unsteady vorticity. In addition, after applying the operators 
— pA 2 V • ( ) and D( )/Dt to the terms in the momentum (5.8) and continuity equations 

(5.9), respectively, and combining the resulting equations, we find that 


-A 2 V 2 p = 0. 


Equations (5.3), (5.11), and (5.12) indicate that first-order entropic, vortical and acoustic 
perturbations of a uniform mean flow are independent modes of unsteady motion. Entropic 
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and vortical disturbances are converted by the mean flow, and therefore have general solu- 

t ions of the form - n c T/ ., t?) 

s = s(rAZ-Vtt) and < = CM,£ ' • ^■ L4 > 

Pressure disturbances are governed by a convected wave equation, and, as a consequence, 
such disturbances exhibit a more complicated behavior than entropic or vortical disturbances. 

The unsteady perturbation can also be expressed in terms of independent entropic. to 
tional velocity. *«, and irrotational velocity, V?, disturbances. Following Goldstein [Gol.6], 
we introduce the velocity splitting 

v = vp + V*<4, ( 5 ' 14) 


where the rotational velocity is 
and pressure are related by 


divergence free, i.e., V x • v R = 0, and the velocity potential 



(5.15) 


It follows from equations (5.8) and (5.9) that 


= 0 or vr = Vfi(r,0,^ V^t) , 
Dt 


(5.16) 


and 


^ Hr A-V*. 

Dt 2 x 


= 0 . 


(5.17 


Thus the rotational velocity, like the entropy and vorticity, is a convected quantity and t e 
velocity potential, like the pressure, is governed by a convected wave equation The .solut 
to this equation can be used to determine the irrotational component, v, - Vo, 

unsteady velocity as well as the unsteady pressure. < . , . , 

For turbomachinerv applications it is usually appropriate to restrict consideration to u - 
steady ^"ons t'hat L harmonic in time and periodic in the circumferential direction. 

Thus, e.g., we can write 


p(x,t) = Re{p(x)exp(iu;f)} = Re 


OO 

£ Pm{r,0eM tm6 ) 


.m=— oo 


exp(iu;f ) 


(5. IS) 


where p(x) is the complex amplitude and tv is the temporal frequency of the unsteady pressure 

Z conve^ed-quantities, say entropy and rotat.onal velocity, have complex amplitudes of 
the form 


and 


s(x) = 


Vfl(x) = 


£ s m {r) exp(imd) 


Lm=— oo 


exp(i/c ? 0 


VR,m( r )exp(im0) 


exp(izc^) 


(5.19) 


(5.20) 


where k ^ = —u>V^ 
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An analytic solution to the convected wave equation (5.12) for the unsteady pressure 
n be determined using the method of separation of variables [TS62, VTM82] leading to a 
solution for the functions p m (r,£) of the form ' ° 

Pm(r,£) — J2n=0 a m,n[Jm(k m ,nr) + Q min Y m (k m<n r)]exp(Xm.nO 

oo (5.21) 

2In=0 a m,n-Em,n(km,n r ) exp(^ mn £) • 

Here a m Q m and k m>n are constants, J m and Y m are Bessel functions of order m of the 
hrst and second kinds, respectively, E(k m , n r) is referred to as a “characteristic E-function”. 

Xm,n = Pm.n + = (1 - M 2 )~ l \lMuA~ 1 If ^( 1 - A / 2 )&2 ^ _ ^2 4-2 1 _ (5 22 ) 


Am ’ n ^ ^ -Y-M J + y (1 — M 2 )k 2 m n - u 2 A~ 2 j . (5.22) 

where M = VJA is the Mach number of the steady background flow. Note that v m n is purely 
imaginary whenever u/A > (1 - AP)^k m , n . Under this condition, the two-components of 
the m, nth pressure pattern propagate unattenuated along spiral paths normal to the lines 
Ktyt + mO+ut = constant. If U /A < (1 - M 2 )k m , n , Xm<n is complex and an m,nth acoustic 
disturbance grows or decays exponentially along the duct depending on whether the - or + 
S1 gn is applicable. The appropriate sign, - or + in (5.22), is determined by the conditions 
imposed on the propagation of an acoustic wave. For example, an acoustic excitation must 
ei er attenuate as it approaches the blade row, or propagate carrying energy toward the 
a e row. us, for a subsonic mean flow (M < 1) the - sign must be selected to describe 

dLnTtTfam eXCltatl0ri ^ UPStream; the + sign ’ to deScribe an acoustic excitation from 

The constants, k m , n and Q m , n , in Equation (5.21) are determined by applying boundarv 
conditions at the inner (hub), r = r H , and outer (casing), r = r D , walls of the duct. The 
so-called ‘hard wall conditions, 


oT = 0 at r — rjj and r = 


(5.23) 


will be used here, but “soft-wall” boundary conditions, in which an acoustic impedence is 
specified along the walls, could also be applied. Conditions (5.23) lead to two simultaneous 
equations for determining the k m , n and Q m , n , i.e., 

Jm(pk m ,n'T'D) + Qm,nY^ {pk m ^r D ) = 0 ( 5.241 


Jm{k m , n r D ) + Qm,nYm{k mtn r D ) - 0 , (5.25) 

where p = r H jr D and the primes indicate differentiation with respect to r. This system has 
a nontrivial solution if and only if the determinant of the coefficients is equal to zero. i.e.. 


^m{p^m,n r D) (/^m,n r £>) 

J'm(km,nr D ) Y^{k m<n r D ) 


(5.26) 
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CHARACTERISTIC NUMBERS 
Ordered by Magnitude 


m n 
1 0 
2 0 

3 0 

4 0 

5 0 
1 1 
2 1 

6 0 

3 1 

7 0 

4 I 

8 0 

5 1 

6 1 

7 1 

Table 5.1: 
increasing 


km.nT'D 

m 

1.353429 

i 

2.682314 

n 

3.957204 

2 

5.178097 

3 

6.335176 

8 

6.576327 

12 

7.060908 

4 

7.462803 

5 

7.844923 

13 

8.560980 

9 

8.830198 

6 

9.639520 

14 

9.987278 

7 

11.232715 

10 

12.507607 

8 


D 

m n 

12.702435 

11 1 

12.826052 

16 0 

12.951395 

9 1 

13.342711 

12 1 

13.802131 

1 1 3 

13.878820 

2 3 

13.896018 

3 3 

14.542588 

10 1 

14.927908 

4 3 

15.067202 

13 1 

15.395325 

18 0 

15.975153 

5 3 

16.277517 

6 3 

16.305273 

11 2 

17.336420 

14 1 


k m ,n r D 

m 

17.518803 

19 

18.063198 

7 

18.444412 

8 

18.689379 

15 

18.907080 

i 20 

19.077499 

12 

19.390278 

9 

19.682482 

16 

19.707966 

13 

19.856277 

10 

20.144039 

17 

20.197456 

14 

20.736031 

11 

20.896011 

1 

21.001701 

2 


Roots of the determinant equation (5.26) with /* - 0.5 
magnitude (only the first 60 are shown). 


Tl 

0 21.182236“ 

3 21.382603 

3 22.127346 

0 22.135458 

0 22.219130 

2 22.223667 

2 22.999720 

1 23.256639 

2 23.469107 

2 23.845089 

0 24.381201 

2 24.771610 

3 24.823000 

4 25.190264 

4 25.360685 

arranged in order of 


, , _ / k ofi\ has a countably infinite number of distinct roots /cm.n^D 

The trancendental equation (5.26) h . in order of magnitude and numbered 

for each integer m, ^whic: can e p ac have n zer os in the interval r„ <r <r D . 

n = 0, 1,2, • • • • The functions E m , n { m ' n > . . • in [VTM82]. It involves a 

A solution scheme for c f the roots, along with a bisection algorithm 

combination of a fixed point iteratio iteration routine or occur in regions in 

for calculating roots .that arc !«t! r 0cdered results of this algorithm are shown 

which the iteratton has drffl^ty m^con ^ ^ ^ constants Qm „ have been 

in Table o.l for // O.o. E-functions, can be determined. A few 

evaluated, the radial mode s “J* 8 t j ue are shown in Figure 3. 

of these functions, normalized by their larg t ic perturbation of a uniform mean 

As indicated by the foregoing results, a gen e d I .const c pe conta ining m 

flow in a cylindrical duct consists of a combination a , mth order Bessel 

lobes and having a radial dependence es the mov e axially away 

paths, unatLuated in amplitude. 


5.2 


Unsteady Perturbations of Nonuniform Mean Flows 


unsl"^ 
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icitv nor to the standard convected wave equation for the unsteady pressure. The entropy 
~ si a purely convected quantity with general solution given by (5.13): but vortical and 
Moust.c disturbances are not independent modes of unsteady motion, and it is not possible 
„ determme analytic solutions for these first-order unsteady flow variables. However bv 
assuming an exponential dependence for the temporal, circumferential, and axial behaviors 
the unsteady perturbation quantities, the first order momentum and continuity equations 

the :r n d * TT but com P llcat «i. second-order, ordinary differential equation for 
he radial behavior of the pressure. This equation along with the boundary conditions at 

he duct walls provides an eigenvalue problem for the unsteady perturbation, which musj 

technique for determining this eigenvalue problem for a mean flow 
with axial shear and swirl is given below. 

VVe consider a steady background flow containing both axial and swirling components of 
velocity but no radial component. The steady velocity and thermodynamic properties are 
described by equations (5.1) and (5.2), respectively. If we assume that the first order urn 
,c“ d} , r ° Py ’ Vel0 T V “ d pressure CMlsist ° f a superposition of fundamental disturbances 
IxtallrecttomeT^ y “ t,me “ d “ circumferential direction and exponentially in the 


= Re{p(r) exp[x(r)£ + i(m9 + ut)}} , 

and set 

A = i(uj + mr'Wg) + , 

the unsteady equations (5.3)— (5.5) reduce to 


As = 0 , 


Aw - 2r-'V,v, + r-'Vfs = -p- (± + ^ + (prA'p'Vfp, 


^ve + ( r l V e + 


dVs \ ■ ,-x-l 

— \v T = —im(pr) V, 


. dV< 

™ ( ' + —v T = -p X p, 


and 


(5.27) 

(5.28) 

(5.29) 

(5.30) 

(5.31) 

(5.32) 


A (pA 2 )- 1 p + (rA 2 )- 1 V g 2 v r + 


d ty 
dr 


+ r V 


dx 

+ ^T^' r + imr ve + \vz = 0. (5.33) 


The unsteady solutions will be composed of convected (A = 0) and nonconvected (A =£ 0) 
components. The entropy is a convected disturbance with axial wave number 


K i ~ -iX = -(w + mr 1 V tf )V^ 1 , (5.34) 

which in this case, varies with radial position. Other convected disturbances can be deter- 
mined as solutions of the momentum and continuity equations. For example, the assumption 
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e *:,i Til and axial (5.32) momentum equations results in 
p un.es, the steady flow satisfies the condition 

(5.35) 


lf . dV e \ / 

-'v # + ^r)(w 


+ mr 


_1 Ve) = -mr 


Vi 


,dV e 


dr 


• ♦ f fled then v = p = 0, and the radial momentum (5.30) and continuity 

If (5.35) is not satished, then v T P _ 

(5.33) equations have solutions of the form 

mVgV^s 


vo = 


V* 


s and — 


2(w + mr 1 Vo) 


(5.36) 


lf , in addition to X = 0, there ate no 

disturbance velocity v< will also be ® umber 0 £ this disturbance is zero. This 

whereupon sy is arbitrary, an e is y mear ly dependent upon r, e.g., if W - *- r - 

condition can be satisfied if the swirl y 1 d P esignat ed as “solid-body swirl,” 

represents the flow behind a high worh 

““If Ze is no axial shear in the mean " ^^il 

is satisfied for the cases of solid-body swirl i a fan rotor (Ker70). Under 

described by Ve = r / r and aPPr ? X1 ^ C u S m f e rential and axial momentum equations become 
either of these circumstances, ... v OT p C an be specified independently, 

iimple multiples of each other, an r , s f ea( fy velocity is zero (Ve = 0) and 

Finally, if the circumferenti a. component of '« ^ 3) that 

m ± 0, it follows from (5.31) and (5.32) that u r a y 


SI 


= m(ojr) 1 V(Vb ■ 


(5.37) 


This solution describes an unsteady velocity disturbance that carries no pressure and is 

convected at the mean flow velocity. q. For such disturbances 

Nonconvected t ' a ^ usly * gCe *, u, and v f in terms of , and 

equations (o.3U) - (O.o-; u 

dp/dr, as long as p _ * +Jr -, V , ( r -. V , + f) * 0 , 

and the results can be subtituted into (5.33) to S wHch 

for p(r). It can be shown that * ^ ^ ^ unstea dy 

all unsteady variables are zero. If D * 0 and we 
velocities are given by 


y r = -\{pD)~ l 


_i f dp 


_ Ver'^VgA' 2 - 2tm(Ar) l ]p 
dr 


vo = {pD) 




Vfl + ~J~ At 


V/A- 2 (r-'V, + ^) + imX 


(5.39) 


(5.40) 
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and 


v t = ( P D ) 


-i^L 

dr 


~^~ r 'Vs {V e A 2 - 2im(Xr)~ 1 ^ p 


~ .\(pA) p . (5.41 ; 


After substituting these expressions into the continuity equation and performing the 
sary algebra we find that 

d2 P , 


+ /i( r )— + / 2 (r)p = 0, 


dr 


neces- 


(5.42) 


where 


and 


Mr) = P^- l Dj^[\{pD)- l }-r- l Ve{V e A- 2 -2im\- l r- 1 ) 

-im(\r)~ x ^r _1 Vg + + r -1 (A“ 2 V 5 2 + 1) - A _1 X^p 


fiix) — — pA A( pDr ) 1 Vg(V$A 2 — 2im\ V *)j 


(5.43) 


-r- 2 V e {A~ 2 V e 2 + 1){V 2 A~ 2 - 2zmA -1 r -1 ) 


-HmA l r 2 



„-i 


V e + 


w 

dr 


+ imX 


(5.44) 


+T " 1Ve lfr ( V * A ~ 2 ~ 2zmA " lr ' 1 ) + A-W - A 2 A- 2 ). 

Equation (5.42) is a rather complicated ordinary differential equation which must be solved 
numerically subject to boundary conditions at the hub (r = r H ) and duct (r = r D ) walls. 
Nontrivial solutions of the resulting eigenvalue problem determine the radial modes p(r) 
corresponding to the eigenvalues x, which determine the attenuation constants and/or the 
axial wavenumbers of the nonconvected disturbances. 

It is interesting to compare the foregoing results to those arising from this same process 
carried out for a uniform steady flow. If V, \ is a constant, Vg = 0 and A = 0. It then follows 
from (5.30) - (5.32) that p(r) = 0, the velocity of this pressure-less disturbance (i.e., the 
rotational velocity) is convected by the mean flow, and this velocity has a solution of the form 
(5.16). If A ± 0, the velocity components are given by Eq. (5.39) - (5.41) with dV^/dr and 
V e set equal to zero, and p{r) is governed by the second-order ordinary differential equation 

S + r ' I * + (i2 - mJr " 2)p = 0 ’ < 5 -«> 

where 

k 2 = A~ 2 [(A 2 - If )x 2 - 2 iuVtx + ^ 2 ] • (5.46) 

Equation (5.45) is Bessel’s equation of order m and has a solution of the form indicated in 
(5.21), and (5.46) yields solutions for \ °f the form (5.22). 
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5.3 Two-Dimensional Flow in a Cylindrical Annulus 

, • • ,• /= o\ throuffh (5.5) we have assumed that the radial component 

In deriving equa ions ^ ’ we also assume that the radial component of 

of the steady velocity is zero [cf. ( • )!• obtain analytic solutions for the 

the unsteady velocity ,s negltgtble, t.e v ° b for the ^ of mean swirl . 
“ unsteady flow — , he d that^^ ^ , ntW varia M s 

itradv ve olay and pressure, within the cylindrical annulus r. - * < r < r. + dr, 


m tne unsieauv vau^uj r — 

are determined by the equations 


Dv 

Dt 


= -pV Xo p 


m = ~ pA ‘ v "‘ 


where x„ = (r„, (, 9). D/Dt is defined in (5.6), v = v ( e t + v»e», and 

d _! d 


Vxo_e4 ^ + r 


(5.47) 

(5.48) 

(5.49) 


In addition, the radial component of the first-order momentum equation (5.4) reduces to 


r- 1 Ve{V e s-2v e ) = P 1 ( r ' A 2v eP~ 


dp 

dr 


(5.50) 


We can unwrap the cylindrical annulus, set e„ - e 0 , V n V e , v n vg, and d/ V 

r-'d/d9, and regard p ag ? r TnTLTiTqul^ and (5.48) govern 

e „ = e, points out horn i h ^ P ^ turb a tio n of a uniform stream flowing at an 

the two-dimensmnd hn d, Ag for the case of uni form mean flow in a 

angle fi = tan (V^J _ _ . ( or _Ce r ) associated with this 

cylindrical duct, the unsteady vor ici y C oressure is governed by a (two- 

two-dimensional flow is convected y e ’ , velocity can be decomposed into 

dimensional) convected wave equation,^ euns^ ^ rotationa l velocity determines 

Th"d“ ^ ^ mean flow, the velocity potential determines 

the pressure and is governed by a convected wave equatrom tangen tial, ,, 

of the first-order unsteady entropy and rotational velocity ar 

5 (X 0 )= £ «m(»*o)eXp[i(lC{,m€ + Kn.m»/)] 


and 


vr(xo) = £ v Rfm (ro)exp[.(K C ,m€ + *,,mi?)] 


(5.52) 


m = — oo 

_i • t i “rircumferential” wave number of the mth disturbance, an 
where K n , m - mr 0 is h ^ number of the mth entropic and rotational 

K^.m — — V u; K ri,m v n) v t 

velocity disturbances. 
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The solution of the two-dimensional convected wave equation tor the complex amplitude 
of the pressure, as determined by the method of separation of variables, is 

CO 

p ( x o ) = J2 Pm{ro)exp( Xm ^ + lKr,,mV) (5.53) 

m = — oo 

where 

Xto ftm + -d-m ”1" COS Cl , 

= (u?V _1 -I- K Vtm sin Q)/(l - Af 2 cos 2 Q) , (5.54) 

<& = (1 — A/ 2 cos 2 ft)' 1 4m - M 2 8 2 m 

and d m is the principal root of d 2 m ; i.e., d m = \d m \ if d 2 m > 0 and d m = i\d m \ if d 2 m < 0. 

If d m ^ 0, each fundamental solution describes two wave-lilce disturbances which, depend- 
ing upon the sign of d 2 m , either grow or decay exponentially in the axial direction {d 2 m > 0) 
or propagate carrying energy away from or toward the blade row {<f m < 0). The condition 
d 2 m = 0 which divides these two types of behavior is referred to as the ‘‘cut off” or acoustic 
resonance condition. In each case, 4 > o, 4 < o or dl = 0, we can write 

Pm (x 0 ) = p m (r 0 )exp[{h m d m + iM 2 6 m cos fi)f + ik Vym tj] , (5.55) 

where the correct value of h m (i.e., +1 or -1) for (P m ± 0, must be determined by whether 
the pressure wave travels in the positive or negative axial directions. 

The two-dimensional, far-field, conditions account for the effects of mean swirl, and for 
this reason they are often applied (in “strips”) in three-dimensional unsteady flow calcula- 
tions. These conditions do not properly account for the radial dependence of the unsteady 
flow variables. Consequently, solutions determined using this quasi-three-dimensional (or 
strip theory ) approach may be of limited use for forced vibration and aeroacoustic response 
studies in which the excitation frequencies of interest are usually high, i.e., u ~ 0(10). 


5.4 Numerical Far-Field Boundary Conditions 


The development of numerical far-field boundary conditions typically involves an exam- 
ination of an approximate set of linear unsteady equations followed by a decomposition of 
the unsteady disturbances into incoming and outgoing waves. The incoming and outgoing 
disturbances at the far-field boundaries are determined by the characteristics of the gov- 
erning equations. The characteristic values or eigenvalues represent the velocities of the 
fundamental disturbance waves, and the eigenfunctions describe the shapes of these waves. 
Conditions must be formulated to neutralize undesired incoming disturbances. Such condi- 
tions are termed “non-reflecting.” 

We can write the first-order unsteady equations (5.3) - (5.5), in matrix form, i.e., 


dn p du p 

sr +A -fr 


^ , du p du p 

TBr' 1 — + C^ + Du p = 0 
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dt 


(5.56) 
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' Vi, 

0 

0 

0 

0 1 

0 

V* 

0 

0 

0 

0 

0 

V* 

0 

r x = 

0 

0 

0 

V* 

0 

0 

0 

pA 2 

0 

Ve J 


'0 0 0 0 0 

0 0 0 0 p~ l 

A= 0 0 0 0 0 , 

0 0 0 0 0 

0 pA 2 0 0 0 

Vi 0 0 0 0 

0 V* 0 0 0 

1 , C = 0 0 Vc o 0 

I 0 0 0 Vz p- 1 

\ 0 0 0 pA 2 V t 


(5.57) 


0 -(pr) l A 2 V 2 


r~ l Vg + dVg/dr 

dVJdr 

pr-'iV? + A 2 ) 


(5.59) 


and r o 0 0 0 0 

r -iy2 0 —2r~ x Vg 0 —(pr)~ x A~ 2 Vg 

D = 0 r~ l V e + dV e /dr 0 0 0 (5-59) 

0 dVJdr 0 0 0 

0 pr-\V e 2 + A 2 ) 0 0 0 J 

The components of the column vector u- are the primitive unsteady flow variables, those of 
the matrices A, B, C and D are determined by the mean or steady backgroun 

One-Dimensional Conditions 

We assume that v r ~ OU 2 ), the unsteady entropy, circumferential velocity and pres- 
sure are related by (5.50), and the circumferential variations (0 / 39 ) m ‘the ^ea _y ^ 
properties are of 0(e 2 ). Then, the components of the column vectors Ad / first-order 
Br-'du? /3Q are at most of 0(t 2 ), and the system of equations (5.56) for the first o 

unsteady flow variables can be approximated by 

^! + C — =o ( 5 - 6 °) 

3t 3i 

Note that, although we have assumed tv ~ 0(e 2 ) in deriving (5,60), a convective equation 

for u r has been retained within this system of equations. 

The matrix C can be diagonalized by employing a similarity transformation, i.e.. 

S _1 CS = A ( 5 - 61 ) 

, <5 rnntams the eigenvectors of C as its columns and A is a diagonal matrix whose 

elements are the eigenvalues of C. Premultiplying (5.60) by S' 1 and transforming to t e 
characteristic variables c = S 'S’ results in 

»£ + A* = 0 ( 5 ' 62 > 

3t 

Equation (5 62) leads to an application of “one-dimensional” boundary conditions to the 
three dimensional unsteady problem. Such an application is effective if the unsteady distur- 
lancet“ves travel along hnes that are nearly parallel to the f-axis. The matrix S , used 
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in transforming from characteristic to physical variables via the equation u p = Sc. In the 
present application S and its inverse, S _1 , are given by 

'1000 0 1 [" 1 0 0 0 0 

01000 0100 0 
S = 0 0 1 0 0 and S _1 =0010 0 

0 0 0 1 1 0 0 0 1/2 (2 pA)~ l 

_ 0 0 0 pA -pA J [ 0 0 0 1/2 —{2pA)~ 1 

The five eigenvalues of the matrix C are 

Aj = A 2 = A 3 = , A 4 = + A , and A .5 — — A . (5.64) 

For subsonic flow, the first four eigenvalues of C are positive and the fifth is negative. 

Therefore, at the inlet boundary, there are four waves entering the computational domain 
which must be specified, and one wave exiting which is extrapolated from the interior. At 
the downstream boundary, waves corresponding to the first four eigenvalues are exiting the 
domain and are therefore extrapolated, while the wave corresponding to A 5 must be specified. 
For example, if it is assumed that there are no external excitations, then the conditions 
Ci = c -2 = c 3 = C 4 = 0 are specified at the upstream boundary, and £5 = 0 is specified at the 
downstream boundary. These conditions imply that 

s = 0, v T = 0 , vo = 0 , 64 = c 5 , and p = -pAc 5 , (5.65) 

at the upstream boundary where the value of c 5 is extrapolated from the interior of the 

solution domain. At the downstream boundary 

s = ci , v r = c 2 , vo = c 3 , vt = c 4 , and p = +pAc 4 , ( 5 . 66 ) 

where, here, the values of Ci,c 2 ,c 3 , and c 4 are extrapolated from the interior. 

Two-Dimensional Conditions 



The one-dimensional boundary conditions are useful if all unsteady disturbance waves 
travel nearly in the positive or negative axial directions. Two-dimensional conditions which 
allow for both axial and circumferential variations in the unsteady disturbances can be 
determined by analyzing the coefficient matrices B and C in (5.56). The approach is de- 
tailed in [Gil90], and its application to three-dimensional flow which yields “quasi-three- 
dimensional” boundary conditions is discussed in [SG91]. 

In developing the two-dimensional far-field conditions, the equations of motion (5.56) are 


approximated by 


du p 1 du p ~0u p 

— + Br - ] — + c— = 0 

at dB 


(5.67) 


This approximation is equivalent to setting v T = 0 in the circumferential and axial momentum 
equations in (5.4) and in the continuity equation (5.5), assuming that (5.50) holds, and 
replacing the radial momentum equation in (5.4) by Dv T /Dt = 0. If we consider a wave-like 
solution of the form 


u p (G 0,t) = u exp[x£ + i(md + u;f)] , 


(5.68) 
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where X ,m and u; are constants, and substitute (5.68) into (5.67), we find that 

[iwl + + X C]u = 0 , ( 5 - 69 ) 

where I is the identity matrix. In general, w and m are specified, and x must be determined. 
A right eigenvector u* of the system (5.69) is determined by solving the equation 


[*wl + imr -1 B + xC]u* = 0 


(5.70) 


indicating that u* corresponds to the eigenvalue -iu of the matrix + xC]- Pre- 

multiplying (5.70) by C" 1 yields 


C _1 [iwl + imr *B + X c l uR = 0 


(5.71 


that u* is also the right eigenvector of the matrix [iwC 1 + imr 1 C *B] corresponding 
to the eigenvalue - X . The perturbation vector u can be expressed as a linear combination 
of the right eigenvectors u* , n = 1.2, . . . ,5, with eigenvalues - X «, >-e„ 

(5.72) 


u p = a „ U n ex P[Xn£] j exp[z(m0 + ut)] . 


Analogously, a left eigenvector u L of the matrix [imr l B + xC] corresponding to eigen 
value -iou is defined by 


(5.73) 

and a left eigenvector v* of the matrix iu,C- + irm-'Cr'B corresponding to the eigenvalue 
is defined by V b c - M + ilBr -*B + *C] = 0 (»•«) 

The right and left eigenvectors corresponding to the eigenvalues X'.- ere thus u n and v„ - 
u^c fcln be shown that the scalar product v‘ ■ »* = 0 for m + », so premult.ply.ng 

n 

(5.72) by results in 

(5.75) 


+ imr 1 B + x^l — 9 •> 


• u p = a n (v„ • uf) exp[xn£] exp[i(m0 + ut)\ , n - 1, 2, . . . , 5 


Nonreflecting boundary conditions are found by identifying which of the 5 
to incoming waves at the boundaries, and spec.fy.ng the constants o„ m (5.<5) tor 

Wa 'por the two-dimensional system (5.69), the eigenvalues are found from the determinant 
of the coefficient matrix, which yields 

[i( w + mr-'Ve ) + Xn^] 3 |[t(w + mr^Vg) + Xn^] + A 2 {m 2 r 2 - X n )} = 0 (5 ’ 76) 

This equation is a dispersion relation that relates the wave number of a given disturbance 
to its frequency. Three of its roots are identical, and are given by [cf. (o.3 )] 

X 1|2 ,3 = — i(ui + mr~ l Vg)Vf l = i*z (5-77) 
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IndThltVhtv t e dist “ rbances “ socia ‘ e J d with the * are convected by the mean flow 
and that they do not decay w,th ax.al d.stance. The remaining two roots are found bv 
solving the quadratic m (5./6), which gives, 


( A 2 ~ ^ 2 )X4,5 = i(uj + mr l V e )V i =F (A 2 - V^)d„ 


where 


(5.78) 


{A 2 - V 2 )d m - AyJiA 2 - K/)(mr>)2 - (w + mr-if/*) 2 . (5.79) 

These results become identical to those given in (5.54) if we replace mr' 1 bv k, and V $ bv 
Vt >- 1 he corresponding right eigenvectors are given by 


if 


' 1 ' 


r 0 ' 


0 


0 

► , uf = < 

1 


0 


0 0 

0 0 

C 

II 

(mV^) _ 1 r 7 

1 

’ > 5 = < 

{ 0 J 


0 j 


. 0 J 



where 7 — u + mr 1 V$. The left eigenvectors v£ are given by 


0 
0 

im(pr )~ 1 > . 

P~ l X 4,5 

-{h + Xa^) J 

(5. SO) 


t l ={1 0 0 0 0 } 


V2 ={01000} 

V3 ={0 0 1 (pVfc)" 1 } 

and 


(5.81) 


'4.5 


= {0 0 im(pr) ' -i-,(pv t )-' WV[[4 ! -t?fe 4J -i7T ( J}. 


The first eigenvector is associated with an entropic disturbance, the second and third with 

ikdTstrbIncei 1 ”” 8 ’ ^ ^ “ d d ° WnStream and upstream moving aeons- 

Assuming that there are no specified external disturbances entering the domain it is 
necessary to specify 


(5.82) 


v n ' = 0, n = 1,2, 3,4 and v£ - tP = 0, n = 5 

at the upstream boundary and downstream boundaries, respectively. 

Alternative Techniques 

Both the one-dimensional and the two-dimensional far field conditions presented above 
prov.de only approximate representations of the three-dimensional, unsteady flow behavior 
m the far field of a blade row operating within a cylindrical duct. Since the full system of 
equations contains variable coefficients, an extension of these techniques to three dimensional 

H , OW , S Z 5 ” 0t readlly a PP^ent. One alternative approach used in flutter calculations is to 
stretch the computational grid in the axial direction so as to dissipate acoustic response waves 
that travel away from the blade row and thereby prevent the occurrence of any reflected 


41 



- . , ran not be used in forced vibration or noise generation studies, m 

sirS-i. -1 '«■“ ”* — — • *• 

-ss — 

able for uniform axial mean flows. Thi mean swirl As suc h, it is only useful 

of acoustic disturbances, but ignores unsteady flow ahead of the first blade row of a 

in special situations, e.g., to re P re ^ , i in modern computational simulations is to use 
machine. The technique genera y P ? flows in a strip-theory fashion. This 

the analytic solutions available for two- ^ the effects of mea n swirl, but it does 

quasi-three-dimensional approa [ 1 unste ady flow quantities. An examination 

31 indicates that this 

approach may not be ’ flows is due to Hall and 

Another approach that has been eouations are first discretized and an eigen- 

Clark [HC93a], wherein the fluid d Y^ d ^ digcretized equations. This method corn- 

value analysis is performed on the resu g ^ tem? and any incoming distur- 

putes the eigenvalues and eigenvec ors o order s of the matrices involved in this 

bances are then specified at intensive. The technique also 

type of analysis are quite large, them M ^ ;t coatains the same difficulties for nonun,- 
requires constant coeffici ’ haracte ristic-based methods described herein, 

form mean flows as the analy ica fi ,. useful three-dimensional solutions for far- 

One potentially promising approach to ® thods for so lving the ordinary differential 

field unsteady disturbances is to in . nonuniform mean flow, 

equation (5.42) for the rad.al behav.o of the mmstead y pr at the duct walls. 

By integrating this equation numen ally subfect to boun * ^ prIncipl e. the 

" SSZ5ZL t numerical solut.on to complete the 

description of the unsteady flow field. 
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6. Numerical Model 


The governing equations for the nonlinear three-dimensional unsteady flow problem were 

Lm7 fs4 Tn this S°T three - dimensional st «“ly linearized unsteady Sow prob- 

„ ” ’ . th “ S«t,o n , approximate representations of these equations will be developed 

g mte difference and finite volume discretizations. Also, the particular finite vol 

ume modei used m the development of the NPHASE and LINFLUX codes will be described 
The numencal model for solving the three-dimensional unsteady, steady and linearized 
unsteady flow problems is based on the nonlinear time-marchinp- FnW 1 • 

developed by Whitfield, Janus and Simpson [WJS88], and 

machinery unsteady aerodynamic applications by Huff, Swafford and Reddy [HSR91] The 

a ed ShTsIT ; n t imPUd ^ mUlti - blOCk ’ finite - TOl “™ * and cimputer code 

cascades This 1 1 P ’ 1 h two ' dime " sion al, unsteady flows through vibrating 

cascades. This analysis has been shown to be robust, with a proven capability to predict 

~ °^ S Wlt sh * r P shock defi mtions, and accurate, allowing second and” third-order 

in rSLH+94Und iT' f deta,Ied . de f cr ;P tion of the latest version of NPHASE can be found 
. , j-n . 6 ref ^ rences Clted therein. Under the present effort the NPHASE code 

been modified to take advantage of the efficiencies inherent in a linearized unsteady 
aerodynamic formulation. The linearized version of NPHASE is called LINFLUX. 

6.1 Discretization of Physical Domain 

In order to develop a discrete approximation to the governing equations the physical 
omain must first be described as a set of discrete points, that form a computational mesh 
or gri . For the NPHASE and LINFLUX analyses, the computational mesh is a sheared 
mes , ypica y generated using the IGB grid generation package of [BH92], This struc- 
tured mesh defines a curvilinear coordinate system, the boundaries of which lie along the 
boundaries of the physical domain, such that there is a one-to-one transformation from 

transf ^ domam ’ ( x )> to P<>mts in the computational domain, (a). Such a 

transformation allows clustering of grid points in regions where flow variables undergo hNh 

Ltion (Ft) “tefW Th ^ P 01 f, m ° tion -. A time-dependent coordinate transfor- 
’ j ’ ' ^ rom the P h - vsical domain, in which the grid deforms with the blade 

motion to the computational domain, in which the grid is fixed, can be applied to simplify 

s u3v detodT h nU T ka ' di ^ erenc ' n S boundary conditions. The numencal grid 
usually defined to be uniform, orthogonal and stationary in computational space for con- 
venience in defining finite-difference approximations. This type of coordinate transformation 

is described by Anderson, et al. [ATP84] and illustrated in Figures 4a and 4c. 

Finite- Difference Coordinate Transformation Metrics 

For finite difference approximations, expressed in the computational domain, we consider 
transformation (X -> A ) I of the governing equations from physical (x,t) to computational 

! : J coordinates - Note that time is defined to be the same in the two coordinate systems 
but has been given separate symbols, * and U, to avoid confusion. For arbitrary phyZi 
omains the transformation from physical to computational coordinates must, in general. 
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be defined numerically, and is typically done by solvmg a set c) This 
equations, as described by Thompson, > ; HM ^^ed to describe the location 

inverse transformation, A-.X, contarns all the^ml and f, Xt/daj , j, k = 

of the grid points in physrcal ^“^"rmation, X - A appear explicitly when 

1.2,3. However, the metrics of the for ^ ^ computat ional coordinates; therefore, the 

the forward 

/•/xnr A\n* svstems, i.e., 


and 


r 

dt ' 


1 

0 

0 

0 

' dt c ' 



dx i 


dx\ 

dt c 

dxi 

5qi 

5x i 
5a 2 

5xi 
5a 3 

da i 


dX = 

dx 2 


5x 2 

dt c 

5X2 

5ai 

5X2 

5a 2 

5l2 

5a 3 

da 2 



dx 3 


5x 3 

dt c 

5x 3 

5ai 

5X3 
5a 2 

5x 3 
5a 3 . 

da 3 

- 

r 

' dt c ' 

r 

1 

0 

0 

0 

* 

df 

■ 


da\ 


da i 
1 dt 

5ai 

5xi 

5ai 

5X2 

5ai 
dx 3 


dx 

1 

dA = 

da 2 

— 

da 2 
'dt 

5a2 

5xi 

5a 2 
5x 2 

5a2 

5x 3 


dx 2 


da 3 

L. 


da 3 
dt 

5a 3 

5xi 

5a 3 
dx 2 

5a 3 
5x 3 . 


dx 3 


= I ax dA 


( 6 . 1 ) 


= J xa dX 


( 6 . 2 ) 


"here J" -d I~ are the matncj o, (the fa ^ Jp. at e inverses of each other, 

lx t) transformations, respective y 

> ’ h . (6.3) 

l.e., jrQ _ (X“ x ) _1 • 

Equation (6.3) provides the ^7/“ ^*r) = (•'“)" 

£ SZSLZZL referred to as metric facob.an. 

Finite Volume Geometry 

For a finite volume discretisation, 

equations, the geometries of the mesh “ .^^xahedral cells which fill the physical solu- 
"aifcrf^itS constant computational coordinate, so that each ce 
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is bounded by the six surfaces, say a x = 7 , 1 - 1 and a 2 = J,J -l and a 3 = K, I< - 1. The 
connectivity of the cells is thus known from the computational coordinates, with neighboring 
cells given by changing a computational coordinate by one. 

The geometric properties of the cell are computed from the locations of the cell vertices. 
For example, the area vector of a cell face is given by the cross product of the directed line 
segments that cross the cell face from “lower left” to “upper right” and from “lower right” to 
“upper left”. Let subscripts correspond to mesh indices, then the area vector of a constant 
a i face is 

A/ = - (x/,y A - - x/, ./_!,*-_!) x (x/,j_i iA - - , (6.4) 

and similar expressions exist for constant a 2 faces, Ajr, and constant a 3 faces, A a -. The 
components of the area vectors can be treated as a matrix, with A jk being the area of 
the constant ol 3 face projected in the x k direction. These area vectors are illustrated for a 
two-dimensional geometry in the accompanying sketch. 



Two-dimensional cell with face area vectors. 

The cell volume can be computed from the six face area vectors and the location of the 
mid-points of the faces. Let the mid-point of a constant cti = 1 face be x Arnid , so that the 
distance between mid-points across the cell in the c*i direction is 


Ax/ x Andd x/_i tin i d (6.5) 

The cell volume is then 

^ = g [(A/ + A/_i) • Ax/ + (Aj -f A d _x) • Axj + (A a - -f A a -_i) ■ Ax a -J (6.6) 

which is equivalent to adding the volumes of the pentahedrons defined by each cell face and 
the center point of the cell. 
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6.2 Nonlinear Analysis 

‘'tn ^of '^s^^o^t.^Fo^c^cad^Tn^hfch^he^Ude^^re^ibratUiglt^e 

row have been given § _ v 4 - 7?(x t) is the instantaneous position of a 

points, x, move in physical space, i.e., Currently two methods are used for de- 

nioving field point having its re erence posi ^ ^ isjpHASE code originally contained 

termining the unsteady i. bJed on an algebraic description. As 

the grid motion formulation of Huff [tin b , d Qn a disp l aC ement field that is 

part of the present effort a grid nw ion or ’ j)i r i c hlet boundary conditions at 

determined as a solution of Laplace s equa ion su NPHASE The latter formula- 

the buae surfaces and in the £ ^SLlation In that larger 

amplitudes S bl^ motion can be considered without grid Une intersections occurnng. 

Finite Difference Approximation 

j * fiTono_rl<ar>pndpnt coordinate transformation, (x, f) * 

th^comp^^rfoned^ordinat^ atd S “lying the chain rnie of part.al different.atron. 

Le " a_a_,a_2iA_ A.- d ^A- . (6-t 

dt 


dt c ' dt da } 

leading to equations of the form 

( d doij d 
\ dt c + dt dcxj 


dx k dx k da 3 


doij d j. a 

u + a^ Ft=s 


(6.8) 


Equation (6.8) can be written in strong conservation form by multiply all terms by 
{J rc t)-i = /« and making use of the tensor identities 


JL ( r r — 

da 3 V dxk > 


dl°* d ax 

= 0 and -gjr+asjp 



= 0 


(6.9) 


After applying the necessary algebra, we find that 


j + dLtj = s 

dt c doij 


(6.10) 


where 


tl = /-u, F, = / c 


^-U + 

k dt dx k , 


and S = / QI S 


( 6 . 11 ) 
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Finite Volume Approximation 


In order to apply a finite volume discretization, the Euler governing equations (3.37) must 
be written in terms of the cell geometry. For example, in three dimensions the cell geometry 
is defined by d, the cell volume, t?j, the volume swept out per unit time by the constant a 
face as the cell interface moves, and A jk , the area of the constant a 3 face projected in the x k 
direction. Thus, we can write the finite volume spatial discretization of equation (3.37) as 

= -8jtj + S = -R (6.12) 

a 

where 

U = tfU , Fj = -tfjU + Ajk. F k , and S = i)S (6.13) 

In equation (6.13) the vectors U and S represent averages over the cell volume of the con- 
served quantities and the source terms, respectively. In equation (6.12) the vector F ,• is the 
flux across a constant a 3 face; the symbol S 3 represents the difference across adjacent cell 

interfaces, i.e. 8 3 ( ) = ( )j+i /2 — ( )j-i/ 2 \ and R is referred to as the residual. The 

repeated j index implies summation over all computational coordinate directions, so that 
the term SjFj is the net flux through the cell. Since the grid deforms as the blades move, 
the cell geometry terms, d, and Ajk, are time dependent. 

The finite difference and finite volume discretizations are seen to be related if the transfor- 
mation metrics of the finite difference discretization are interpreted in terms of cell geometry. 
For example, in three dimensions d is analogous to I ax , i) 3 is analogous to -I ax da 3 /dt. and 
A jk is analogous to I ax dcx 3 jdx k . This interpretation shows the similarity between the finite 
difference definitions in equation (6.11) and the finite volume definitions in equation (6.13). 

In the NPHASE analysis [HSR91], the fluxes at the cell interfaces are computed from 
the values of the state vector on either side of the interface using the approximate one- 
dimensional Riemann solver developed by Roe [Roe81]. The interfacial flux is computed by 
multiplying the difference in the state vector across the interface by a flux Jacobian matrix 
representing the local interface conditions. The eigenvalues of this flux Jacobian matrix are 
used to determine which characteristic modes are included, thus controlling the direction 
of spatial differencing. This technique is known as flux difference splitting and results in a 
first-order accurate spatial differencing scheme. Higher order accuracy is obtained by using 
a corrective flux, which is limited by a TVD scheme [OC84] to control dispersive errors. This 
flux differencing scheme is described in §6.4 for the steady and linearized unsteady fluxes. 

The NPHASE time differencing is based on the three level Beam- Warming [WB77] rep- 
resentation of the time derivative. Thus, if we let the superscript n refer to the current or 
nth time level, we can write 


SU 

dt c 


~ n+1 - n ~ n 

U -U = AU = 


1 + 'P 


0Ai « + 


dt 


dt c 


~ n— 1 

+ 'PAU 


(6.14) 


where 0 and »P are constants that define the relative weighting of the time levels. For 
example, setting 0 = 1 and 'P = 1/2 results in the second-order accurate, implicit, three 
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point backward difference s c h e me currently *the 

“ («•«)• “S -ms and substituting for 0 

and ^ results in _ ^ -i (6 . 15) 


2Af c 

"rttttr^orlV^ceil volume results in 

* / jin-1 \ 

= 0 . (6.16) 


= n+1 

= -R 


f Tn + 1 


-U n + 


2 A t c = n +! 

- 

3i? n+l 


d 


171 — 1 


3i? n+1 


(U n - U” l ) + 


u n 

^n+1 


r 4 tT' 1 ’ 

r+1 _.r + - T 


•6 n+1 _ R n+1 (U n+1 ) is a nonlinear function of U” +1 . 
Equation (6.16) is nonlinear, since procedure in which the residual equation 

Therefore, it is solved using a Newto flux terms are approximately factored to 

is linearized using flux vector sp i mg block- tri diagonal, system of equa- 

ELSE m « - - — the — 

Ration (.6.12, must be ^ £ 

( :, 2S ) at the blade the prescription of externa, urn 

outflow boundaries. The tar- held con disturbance waves coming from within 

steady aerodynamic excitations, and aU» and outflow boundaries without reflection, 

the solution domain to pass rou S . . , iqpHASE analysis was based on characteristic 

The flow tangency condition in the origin imple y mentation . This condition was 

theory, and one phantom cell was us svmme try condition, which lowers the spurious 

changed to a two J at ed at k blade surface. The current NPHASE, 

numerical entropy and vorticity that .s gener ^ ^ previous time step when computing 

multi-block, nonlinear analysis uses th periodic block boundaries. In 

across a block interface, thus s Ltched in the axial direction 

the original NPHASE ana ysis e c current effort, far-field conditions based on one- 

to dissipate outgoing waves. As part and on appro ximate two-dimensional 

dimensional characteristic theory, cf. equatic ( ;60)^ PP for the prescri ption 

characteristic theory [Gil90], were implement^. Th * £ ^ disturbances , if the 

of incoming unsteady excitations, but ^ ioo As a result we used the 

-ter have calculations reported in §7. 

stretched mesh capability tor trie nr an. 

6.3 The Steady and Linearized Unsteady Analyses 

The field equations that govern the ^ v ^ty o' 

flows through a blade row rotating at const § have been written in strong 

as equations (4.46) and (4.47), r«pectively. Th^equat^ upon „ and/or R 

conservation form in which terms that depen P “ flow VMiables u and the displacement 
SdT“ed e toTof T ^) 2“ is a pseudo- time variable. The displacement held 
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must be prescribed throughout the solution domain i e a sinrie evt.n^ ^ u, , 
region of finite extent in the axial direction In addit L 'lh ♦ !’ 6Xt ^ ded ’ bl ade-passage 

field equations must be solved subject to annrnn/ !?. > ^ “ d llnearized unsteady 

due, walls, the blade-to-blade V 

of the computational domain. d outflow boundaries 

Deforming Grids and the Dual Coordinate Transformation 

The equations that describe the steady and linearized unsteadv «™. k , 
in terms of the independent variables /and t where xTs ht f bee " 

position of a moving field noint Thi. d V * ,s ' h referenre or steady state 

™zz, 

points, i.e. x(x t) = v + T?(£ f\A / ! T , he time ' de Pendent locations of the grid 

Figure 4a) which, if R = 0, dbfu/e/forT' ^ 

spa^r^:t«r: : h iz met r ter r ™ to **“ «■- 

time. This dual coordinate transformation, (x,f) _ ( S “.“independent of 

equations for the steady and linearized unsteady flowL th“t ie b^d - -V” 

~/uat,r Ce tomS ’ reSUUinS fr ° m ^ *“ -*-*>« ^ thel“L^ 

For Unsteady' flow/t/oug/wbratfng^l^^rows^^v^tor^^sho^ld 1 ^ 115 /^ Pt escr ihed. 
L h z h X:i:zt a t h“«Te s r*::: dm/ e t: b aini e - * r ^ 

whatever manner is most convenient for implementing the flow bound/ condlns. “ 
The blade motion and, therefore the deform a t inn fijaM u 

instantaneous position , oi a field ^ 

X = x + i?e{R(x) exp(iwi)} . (g. 17 ) 


R(x) = R Bn (x), 


( 6 . 18 ) 
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and. since the duct walls are rigid, we set 


R(x) = 0, at f = r H ,r D 


(6.19) 


.x, , convenience in impos.ng the inflow and outflow boundary conditions, the deformation 
field is assumed to vanish at the far field boundar.es, i.e., 


R(x) = 0 for 




( 6 . 20 ) 


where f. and fr are the axia, — 

, boundary conditions, the deformat.on field R(x, , 

assumed to satisfy the condition 

R(f, e + ‘2ir I N , o S Rfr 0 exp(zna) 

, , , • , , • i|.p interior of the computational domain can be 

The displacement or deformation Laplace’s equation to obtain the smoothest 

chosen arbitrarily, and is defined here y g ^ maximurn displacements on the bound- 

Jd **> we require a solution of the equation^ 

V|R = 0 (6 ’ 22) 

over a single, extended, blade-passage region, subject to the Dirichlet boundary conditions 
(6. IS)— (6.20) and the phase-lagged periodicity condition (6. )■ 

Finite Difference Approximation 

To a. 

ordinates we employ a time-in e P ’ g Because the unsteady grid deformation 

reference physical space to comp transformation, there is no unsteadiness in 

£ =1" r^^afus by applying the chain ru.e for the 

stationary x -*• a transformation, i.e., 

d__ _ do; 0_ (6.23) 

dx k ~ dxkdotj ’ 

a „d man.pmaf.ng fhe resuming expressions info strong conservation form by us.ng fhe tensor 

identify a_ _ n (6.24) 

3 '.j V dxk) 


da. 


zeroth-order or mean flow ^nation 

becomes 


dv 

Or 


d F, 


la + da ’ 


S , 


(6.25) 


50 


where 


U = /"U , F, = I°‘pL F t and S = /«*S . 


The perturbation equation is 


where 


and 


< 9 u 

dr 


, • - , dfj ds _ 

+ zo;u + — — — U = 


df f 


a 


da, du da, 

+ r[-iw(Vx • R)U + (V* • R)S + (R • V*)S] , 


u = I ai U , 


f rax da, r T axdaj 

h ~ 1 oT-Iife = i -zz — ^77U , 
dxk dxk 


ed _ rax da, 
J dx k 


f) R 

—iu>R Xk U — k F m + (Vx * R)F^ 
ox m 


(6.26) 


(6.27; 


(6.28) 


The right hand side of equation (6.27) contains source terms, due to the grid deformation, 
which depend upon the specified deformation field and the mean flow solution. 


Finite Volume Approximation 


In order to apply a finite volume discretization, the integral forms of the governing 
equations for the steady (4.48) and the linearized unsteady (4.49) flows must be written 
in terms of the cell geometry. For example, in three dimensions the mean cell geometry 
is defined by t?, the mean cell volume and A, k , the mean area of the constant a, cell face 

projected in the x k direction. Thus, we can write the finite volume spatial discretization of 
equation (4.48) as 


5U 

dr 


= -6, F, + S = -R 

a 


(6.29) 


where 


U - i?U , F, = Aj k F k , and S = t?S (6.30) 

Here, the vectors U and S represent averages of the state and source-term vectors, respec- 

tively, over the mean cell volume and R is the residual, which should not be confused with 
the complex amplitude (R) of the grid-point displacement field. The vector F ; is the flux 
across a constant a, face, and the symbol 8, represents the difference across adjacent cell 

interfaces, i.e. 8,( ) = ( ) j+ 1/2 — ( ) J _ 1 y 2 - The repeated j index implies summation over 

all computational coordinate directions, so that the term 8,F, is the net flux through the 
cell, w’hich is balanced by the source term once the pseudo-time integration converges. 

In addition to the metrics d and A, k , that describe the mean cell geometry, the pertur- 
bation equation contains terms that depend upon the displacement R, which must also be 
interpreted in terms of cell geometry. For example, in three dimensions Ad ~ d — d is the 
instantaneous cell dilatation, d, = iuAj k R Xk is the volume swept out by the constant a, face 
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per unit time, and the fljk ~ A.jk 
Thus we can write (4.49) as 


— A 3 k are the perturbations in 


the cell-face area vectors. 


— + iu;u + Sjfj 

dT \ a 


J®. u = -tw( A t?)U - Sjif + (A#)S + d{R ■ Vx)S 


(6.31) 


where 


- ~ le a / Ffe - u and i d = -i?jU + . 

u = i?u , f j = u ,ana i, J 


Written in terms of a residual, equation (6.31) becomes 


6>u 

dr 


e th*. left hand side of (6.31) is computed using flux difference splitting. 
The flux term, fj, on the left h \ grid-deformation, source terms, and 

The terms on the right hand s.de o the deformation field. R 

the steady and linearized unsteady equa- 

tions are 


- 8S _ 

= -iwu - 6jfj + 


tw( Atf)U - Sjif + (Atf )S + fl(R ' Vx)S - 


= — r 


(6.32) 


(6.33) 


AU n ,U n+1 - 
A r Ar 


U r 


_R n+1 and 




Au 

Ar 


= _ f „ +1 

Ar 


(6.34) 


respectively, where the mean U^Swton iterations are 

equation is linearized using flux vector sph ^."“bloTk rSw equations. Symmetric Gauss- 

at the blade surfaces and duct walls 1 olade-to^lade P t he unsteady field equation (6.31) 

conditions at the inflow and outflow (4.33), the phase-lagged 

must be solved subject to the flow ^ allow the prescription of unsteady 

periodicity conditio ns (4.8)^ disturbance waves coming from within the solution 

domataThrough the inflow and outflow hcund^ries *‘‘ h °" t is base d on charac- 

The flow tangency condrtron, now t used m ^ “ y ^mentation, which lowers 

described in §5.3 and 5.4. 
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6.4 Discretization Scheme 


In the nonlinear unsteady analyses of [WJSS8] and [HSR91] a cell-centered finite-volume 
method is used to spatially discretize the governing equations. The residuals due to net 
fluxes at the cell interfaces are obtained by flux- difference splitting, obtaining an approxi- 
mate Riemann solution with Roe averaged variables [RoeSl], A corrective flux is used to 
obtain higher order spatial accuracy. This particular numerical approach offers the abilitv 

osc^llatimis dlSCOntinUltieS ^ P ° intS without Producing spurious non-physical 


Evaluation of Steady and Linearized Unsteady Flux Terms 

The flux vectors associated with changes in the steady and unsteady state variables, i.e.. 
j and f, in equations (6.29) and (6.31), respectively, are evaluated using a flux splitting 
procedure, which is based on fluxes normal to the cell faces. The inviscid flux vectors Ft 
are homogeneous in U of degree one, i.e., F k = (dF k /dU)U. Thus, in computational space! 

he flux vectors for the steady and first-harmonic unsteady flows [cf. (6.30) and (6 v>9)l can 
be written as v n 


and 


p. _ 3 *• _ j dfj 

’ ,k “ - Aik W v - du 


U 


f - 4 f _ 4 ^F, dF } 

fj “ A ^ k - du u = du u ’ 


(6.35) 

(6.36) 


so the functional relations F,(U) and f,-(u) are identical. Equations (6.35) and (6 36) relate 
the transformed flux vectors F, and f) to the cell face geometries, A jk , and the dependent- 
variable state vectors, U and u, in reference physical space. Both the steady and unsteady 

uxes are computed by using flux difference splitting to solve an approximate Riemann 
problem, as outlined below. 

The flux splitting is based on a similarity transformation of the computational flux Ja- 
cobian matrix dF ; /dU. We set 

dF 

~q\J ~ (T^IT x )j (6.37) 


where T is the matrix of right eigenvectors, A is the diagonal matrix of eigenvalues, and T -1 
is the matrix of left eigenvectors, of the flux Jacobian matrix, dfydU. From equation (6.35) 
it can be seen that this transformation is determined by the A jk and the physical flux 

Jacobian matrix (dF k /dU). The flux Jacobian matrix is split into the sum of right and left 
traveling modes, i.e., 

dF* 

W = (TA ± T -1 )j (6 . 3S) 


where A + and A are diagonal matrices containing the positive ( + ) and negative (-) eigen- 
values, respectively. Positive eigenvalues correspond to right traveling disturbances °and 
negative eigenvalues correspond to left traveling disturbances. The sign of the wave speed 
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(U „ A+ or A -) determines the direction in which spatial differencing is applied. Thus, the 
flux vectors can split according to 


F .A 

dV 



(6.39) 


r the negative coordinate direction and the ( ) 

with the (+) terms using in °™ a j coordinate direction. 

terms using information from the ion ires that the fluxes at the cell surfaces 
A cell-centered finite-volume dis t <1 ^ ^ volum es. One method of 

be computed in terms of val, res . of the * t Riem ann blem . Ro e [RoeSl] solved 

determining these fluxes is to solve an appmx^ nUng average interface conations. 

this problem by defining a flux Ja “* oooos i te 'sides of a cell interface be denoted by the 
Let the values of the state vector on opposit :ej^ ^ ^ in t er face by the subscript / An 
subscripts L and R, and the va uej^ ^ is defined using the following relations 

intermediate state vector, Hoe, JJIEt l+ 

The flux at the interface is constructed from p‘^tet^‘^«r across 

plus the flux due to waves approaching the interface 

the interface. vectors, F and f, based on U = U L and the 

dac^^: w 561 


F, = F(Ul) + 


aF 

au 


IU = U R oe/ 


(Ur — Ul) , 


(6.41) 


and 


ii = f(u L ) + 


(ur - u l) , 


(6.42) 


y lu=UR 0 e/ 

where the disturbances with negative wave ^speed have the 
rrii "or aJLage of the disturbances traveling at posit, ve 

and negative wave speeds. vec tor splitting proposed by Steger 

An alternative method of flux splitting » the «ux J ^ [actoriza tion of the time 

according to the sign of the wave speed, ,.e„ 

\ + 


Fr = 


I aF 




Ul + 



Ur 


(6.43) 


IU=Uh; 
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Spatial Differencing and Corrective Fluxes 


Once the fluxes have been computed, they are spatially differenced to compute the net 
flux through the control volume. The difference expression is 


= E 

j 



(6.44) 


where j corresponds to the computational coordinate direction, and fractional grid indices in- 
dicate cell interfaces. This difference approximation is first-order accurate if the flux is based 
only on information from adjacent cells, but higher order spatial accuracy can be achieved 
by adding a corrective flux. The corrective flux brings in information from additional neigh- 
boring cells, but requires the use of a flux limiter to control dispersive errors. Various flux 
limiting schemes exist, and the choice between them is not clear. The NPHASE code cur- 
rently supports min-mod [OC84], Superbee [Roe85], and Van Leer flux limiters [VL74], with 
Van Leer being preferred. 

In the present implementation flux limiters are used in the steady analysis, but not in 
the linearized unsteady analysis. The limiters used in the steady analysis are activated by 
changes of sign in the jumps in characteristic variables at adjacent interfaces, such as occur 
at shocks and at stagnation points. The corrected steady flux is given by 


F ;Li = ®i( u i) + I Iff 




dU 


u=u 


Roe 


(U f+i - Uj) + ^T[<t + - <r ] , (6.45) 


where the limited jumps in the characteristic variables are 


" + = L(<r+_i,<r+ + J and tr =L (<t. + l, ov + |)] , (6.46) 

and the jumps in the characteristic variables at the j - f cell interface are given by 

a t-f = yl±T_1 ( U j-f+i “ (6.47) 

The matrices T and A. in (6.47) are computed based on Roe averaged variables, and the 
function L in (6.46) is a generic limiting function. For example, the Van Leer limiter can be 
expressed as 

V >1, ' ■ »' , (6.48) 

C7i -f U 2 

where the superscript n refers to the nth component of the column vector. 

No flux limiter is used for the perturbation fluxes, so the corrective flux is comprised 
of left traveling waves at the upstream interface of the adjacent upstream cell and right 
traveling waves at the downstream interface of the adjacent downstream cell. These waves 
at adjacent cell interfaces are approximated using the Roe averaged matrix at the current 
interface. Thus, the corrected perturbation flux is given by 


$ 






u=u Roe/ 


(Uj+I - U j) 


1 

+ 2 



(6.49) 


u=u Roe , 


(Uj -Uj-i) - 



u=u Roe , 


K+2 - U J + 1 ) 
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This expression results in second order spatial accuracy. 

Evaluation of Grid Deformation Source Terms 

in the linearised unsteady analysis 

the governing equations are expressed m .ret P . locat i on . The volume source 

ar= Sue to changes in the first order, in R, perturbation 
term is given by -zu;(Atf)U + ( ) ’ The cell face-area, source term, f , , accounts 

in cell volume, as determined by .qua»on (6^ mcel ^ ^ ^ (6 . 32) . The 

for the mean flux through the movi g ^ computed using first order expansions m R 
perturbations in projected face a ’ J ’ eauatio n (6.4). The swept volume is given by 
for the area of a cell face ’ as e ^ R is ta ^ e n to be the average over a cell face. T e 

i)j = i*A ]k R: k , where the displacem R 3g)] accounts for changes in radial location, 

remaining grid detonn-tK® ^source l U R based on the average d.splacement 

In finite volume form this term is tffK. 
of the cell vertices. 

Pseudo-Time Integration 

Solutions to the nonlinear steady “^“^^J^^Mlutlons are determined^ 

by integrating in pseudo trme unt.l a com«g ^ impUdt> appro ximately-factored 

The pseudo-time integrations are Perfo « d<rf , me discretization results ... a 

method [VW93], For the steady analys ns, _P ^ ^ both ^ steady and linearized 

nonlinear equation that is solve us ® . ation error f s reduced using Gauss- Seidel su 

unsteady analyses, the approxima ' time d ; scr etization requires evaluations of ux 
iterations, as explained below, T P „ based on Steger-Warming flux-vector 

Jacobian matrices. These are calculated numerical calculation of these 

splitting, as opposed to a Lfuing. A simplification in the linearized 

matrices, based on Roe-averaged flux- ^^^ ^ reference physic al coordinates is that 

^ even though the computati m 
def ^^u^SiS; scheme chosen for the nonlinear steady analysis, cf. (6-34, 

is expressed as L(U n+1 ) = (U n+1 - U n )/0 + R n+1 = 0 * 

r - n TT" +1 since R n+1 is a nonlinear function of U n+1 , and 
where L is nonlinear in U , since rv 


0 = 7 ^ 


At0 

(i + Wfi 


At 

T 


(6.51) 


iteration index, then the iteration formula is 

(6.52) 



\ ( U P - U p_l ) = -L(U P_1 ) , 
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where U° = U n and AU P = U p - U p_1 is the Newton update to the state vector. Once the 
Newton iteration converges, AU P = 0 and U p = U n+1 . The flux terms on the left hand side 
of equation (6.52) are computed using flux vector splitting, for convenience in defining the 
approximate factorization, while the flux terms on the right hand side are evaluated with 
flux difference splitting, corrected for higher order spatial accuracy. 

The Newton iteration matrix in (6.52) is given by 


du + du ’ 


(6.53) 


where the residual, R, is defined in (6.29). The change in the residual due to the Newton 
update thus depends upon a change in the net flux and a change in the source term. The 
change in the net flux due to the Newton update is 




dt p ~ l 

dV 


AU P , 


(6.54) 


and is evaluated using flux vector splitting. To apply flux vector splitting, the geometric 
terms defining <9F/dU must be evaluated at the proper location. For example, consider the 
flux crossing the j - 1/2 and j + 1/2 cell faces, which are the boundaries of the cell in the 
ocj computational coordinate direction. When evaluating flux Jacobians using flux vector 
splitting, let the first subscript indicate the cell index and the second subscript indicate the 
interface index. The flux in the a, direction as evaluated using (6.54) is 


-l 


dt] 

8U 


\j,j+l/2j 


( au ?) + ( tT 


ip- 1 1 


dt 

~d\T 




dU 


w + ( t 


dU 


j+ lj+l/2; 


j-lj-1/2, 


( AU '«) 


K.) 


(6.55) 


which is based on the flux vector splitting given in (6.43). The Newton iteration formula is 
thus 


0 AU P + 6j 


'dF p 

dV 


-i 


-AU P 


^S p 

mT 


-i 


-au p = -e-^up -1 - u n ) - R p_1 


(6.56) 


where the repeated j index for the flux term implies summation over all coordinate directions. 
Once the Newton iteration converges, Equation (6.50) is satisfied. 

For the linearized unsteady analysis, the particular time discretization scheme chosen, 
cf. (6.34), can be expressed in the form 


L(u n+1 ) = 0 ‘(u n+1 - u n ) + r" +1 = 0 _1 Au + r n+1 = 0 (6.57) 


where Au is the pseudo-time update, and L is linear, since r n+1 is a linear function of u n+1 . 
Equation (6.57) is solved by using a Taylor series expansion in pseudo time for the residual, 
i.e., 


p”+i 


= r n + 


dr / 
fa ( U 


n+ 1 



(6.58) 


57 



This can be combined with (6.57) to obtain 

(s-'i + |D) (“ n+I - u ") = - f ” • (6 ' 59) 

K„,h the left- and right-hand-sides of equation (6.59) are computed using 

once, and it .eads to much better convergence rates for 

linear unsteady solutions. change in the residual due to the 

volume form, the net dux term, and 

the source term, i e., 


dr 

du 


Au = iu>i?Au + 8j 


d¥j 

dU 


Au n 


_ , 

dv 


(6.60) 


where the dux term is evaluated using dux d.derence splitting. The pseudo-time update 
formula is thus 


0 _1 Au n + iuh?Au n + 8j 


d l± 

d U 


Au” 


_ ^t-Au n = -r n 

dU 


(6.61) 


where the repeated j index 

steady problem is linear, Newton 1 era 1 .. . which are incorporated into a Newton- 

- be •— in a 

“^ coupling between adjacent 

and (6.61) expensive to evaluate, so e order to re duce the error introduced by 

product of a positive and negative opera [WT9 i] are used as part of the 

tim^st^pi^ a pro^^^ a ^^' approximate factorization equation for the Newton iteration 

procedure used in the nonlinear steady analysis is 

(e-‘i + i) ^ - d.au^au Wi au^ = -e-Hu- -u«)-R^ 

and the approximate factorization equation used for the pseudo-time stepping in the lin- 
earized unsteady analysis is 

(6.63) 


0 -1 I + 


?L) Au n » DjAu” - M+. 1 AuJ_ 1 + Mj +a Au" +1 - 


= -f n 


du 


In bosh equations the j j subscript is index ““eNewmn ite^ion iX.^Thl 

— ZZS2Z - & “ - on the state vector 
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U /~‘ at ‘he previous Newton iteration. The D matrix represents the diagonal elements 
o he iteration matrix. The M+ and M matrices represent the off-diagonal elements 

ot the iteration matrix in the negative and positive computational coordinate directions 
respectively. 

Introducing l as the Gauss-Seidel iteration index the sub-iteration formula is given by 


DAV/ - + MJ+.AU#' = -L-' 

DjAUf +1 + M- +1 AU^f -Mlh.AUjd, = -Id- 


le. 64) 


/^UP./+2 _ — AU P,/ 

The first sub-iteration is over positive grid indices and the second sub-iteration is over neg- 
ative grid indices. The sub-iteration procedure is thus an LU decomposition of the Newton 
iteration or pseudo-time matrix, with forward and backward substitution. Once the Gauss- 
Seidel sub-iteration converges, equation (6.56) or (6.61) is satisfied. 

Once the solutions converge to a steady state, any error introduced bv the Newton 
iteration or the approximate factorization vanishes. Only the error in the residual calculation 
° eqUat * 0n ^-“9) or equation (6.33) remains. These residuals are calculated using flux 
difference splitting with a corrective flux, as given by (6.45) and (6.49), to obtain second 
order spatial accuracy. 
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7. Numerical Results 


aerody T ic e «T n < r " 

: h ; b -:?L a ** ^ 

" cits ^ by ^ .** .»*« 

for the Tenth Standard Cascade operating at an inlet Mach number M and fljjrgle, 

f> of 0 7 and 55 deg respectively, and at an inlet Mach num er an ow a S 

58 Xa In the first cie the steady background flow is entirely subsume; m the setxmd U 

is transonic with a norma] shock emanating from the suet, on surf, tee of F ° r 
, , v flnw<? exc ited bv external aerodynamic excitations, the Tenth Standard Cascade 
“s tp^ing at SI* - o', and ^ = 55 deg, and the steady background flow ,s ent.rely 

SUb r° n ^iaition to the LINFLUX results for the Tenth Standard Cascade, for purposes of 
comparison, we will also present nonlinear unsteady 

the flows excited * as 

linear response PredKtwns for Ud b utter, the unsteady flow is 

d “ ' t rimall pertSl off^alUeadyUground flow. The steady version 

SPHtSE analysis is used to provide the steady background flow information for the 

UNFLUX^ lalcul^ ^ the M-potLial analysis CASPOF [Cas83] is used to prov.de th.s 
LINFLUX calculations tne appropriate unste ady far field cond.t.ons 

LtrnotTe. been iniorporated, at present, the NPHASE analysis cannot be applied to 
predict unsteady flows excited by external aerodynamic excitations. 

7.1 Flow Configuration and Computational Meshes 

The Tenth Standard Cascade has a stagger angle, 0, of 45 deg and a gap/chord ratio G 
o, unity (see F.gure 2). The blades are constructed 
of a modified NACA four-digit series airfoil on a circular 

distribution is given by 

T(x) = Hi[2 Mx'I> -1.26x-3.516x ! + 2.843x=- 1.036x1, »<x<1. (7.1) 

in tedge-shaped trailing edges. The camber distribut.on ts g,ven by 

C(x) = H c -R+[R 2 - (* -0-5) 2 ] 1/2 7 0 < x < 1 , 


(7.2) 


it M 01 is the height of the camberline at mid-chord and R - {21 Ic) {H c + 0.2o) 
r s tie r!d C ius of’the circuliarc camber line. Thus, the surface coordinates for the reference 

(n = 0) blade are given by 

[l , !( |e. ± = [x + 0.5T(x)sm«, C(x)±0.5T(x)cos«), 0<x<l, U3) 
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where 6 = tan -'(dC/dx). The Tenth Standard Cascade is defined by setting Hr = n 06 „ A 

%t A ' ai j d ttlerefore > a ca5ca de of modified NACA 5506 airfoils 

FI ^ ,• h " earlz ' d unstead y solutions were determined, using NPHASE and LIN 

FLUX respect, vely, on H-type meshes. The meshes used in the blade vibration L . 

r 75 points ™ pi ^ ten* 

downstream from the blade row to dissipate X“es ^ ^ “ d 

the 1 CASPOFcJc 1 

extended one axial chord upstream ^d Iwn t^alto^tt wjT °Ti ** ^ 
flows, radial mesh lines were concentrated n^Ss ^ ^ ** *““■* 

370 WorkstaUo^Converge^^PffAS^stea^v^ 111 ' Cal , U !f ltions performed on ™ IBM- 
pseudo-time iterations. Ihe’ steady suCc^luZs « "a 'm Tu \ ^ ^j' 000 

75 P CPU^tate C ThTmxt; tr 1 a, ^ i0n,C solutlons (2« x q 81 mesh), 
excited by blade vibrations or^St"^ 

vibra^ionsj'also performed oTa^S^x ^Tmefh ^ ^ tranS ° niC fl ° WS eXdted by bkde 
The calculations for unsteady subsonic flo^ r C , PU mimiteS ’ 

were performed on a 241 by 81 mesh, also required about 20 CPU m^n 7 ^ WhlCh 
The NPHASE nonlinear unsteady calculaSs were ^st^ frorT Z "*7^ 

solution, and performed using 1,000 time-steps per cvcle of hi rl a PP ro Pnate steady 

The subsonic unused unsteady calculations were performed on a 121 x „ u a 
required about 40 CPU min per blade oas^^ ^ , U1 X 41 H - mes h and 

solutions, performed on a 241 x 81 HmeJ ^ COnver S ed transoni c inviscid 
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j j steady W 1 potential solutions were obtained using CASPOF, ^within 
Tc?V secondhand LINFLO, linearized, inviscid solutions requ.red about 90 
per unsteady case. 

7 2 Blade Vibration 

... astsar 

memo, blade motions at small * rotations (torsion ) about axes at the blade 

normal to the blade chords ( en g) , nant if all ac0 ustic response disturbances a 
midchords. Blade motions are ^nn ^ {rQm thg bl ’ ade row; superresonant (m, n), i m an 
tenuate with increasing axia ,, DS tream and far downstream regions, respective j, 

n such disturbances persrst m the far upst ^ , f ^ ^ one acoustic response 

ituTb're^SstsTn eZ the far upstream or far downstream regions and carr.es energ. 

along the blade row [Ver89b]. -0.7 and ^ = 55 deg and exit 

For the subsonic cascade tential calculation, of M_c* = 0-446 and tt-oo - 

flow conditions, as determine y , degrees) are a_ x = -26.93a;, a ^ 

40.2 deg, the resonant interblade ^ gubscripts refer to the far upstream (-00 

117.12a;, = - 31 . 80 w and <r +0O ■ g Scripts indicate that there are two resonan 

or far downstream (+ 00 ) regions regions. For the transonic configuration 

interblade phase angles associated with each ° f = 40 . 3 deg, the resonant phase 
with M_oo = 0.8, 1U. = 58 deg, M, = 0.432 andJ4^ ^ ^ = 06 j n both 

“ gleS ,h r e e bl^e molio^at’ fnterblade phale £*• lying between the lowest (»„) « 
highest (oij resonant interblade phase angles are superresona . 


Subsonic Flow (M- <*> - °- 7 , B -00 55 deg) 


t , k1 a de- surface Mach number distributions for the 

The steady Mach number con ours an ^ M determined by the nonlinear Euler 

subsonic operating condition ( M -°° . ’ h face Mac h number distributions predicted 

analysis are shown in Figure pred icts a mean flow exit Mach number, 

using the CASPOF analysis. The ^ respectively , and introduces spurious 

M+ 00 , and flow angle, fl+oo, of °. d 0 01 % on the suction and pressure surfaces, 

— * — 9ow ex,t Mach nu 

and flow angle of 0.446 and 40.2 _ r ®*? e ^ bso ^ c NA CA 5506 cascade are shown m Fig- 

Unsteadv response predictions g , ith blade displacement) component of 

ures 7 through 10. Contours of the ou^-ph«ejw> ^ ^ UNFLUX rmalyses. 

the unsteady pressure, . w de vibration at unit frequency are shown m 

“ be in - good agreemen 

tor 

shownln figure ^mSt! of The n^Lear Euler analysis NPH ASE. the linearized Euler 
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analyse INFLUX, and the linearized potential analysis LINFLO are seen to be in excellent 
agreement tor the superresonant (1,1), in-phase vibration, and in excellent agreement all™ 

atnv7h SU M 7 faCe f ° r the ^resonant, out-of-phase vibration. However, differences exist 
. 6 be bbl<ie suction surface between the results of the two Euler analyses and those of 

the potential analysis, for the out-of-phase motion. f 

The blade motions considered in Figure 8 are both stable, i.e., W c < 0 The pressure 
displacement function curves for the in-phase torsional motion reveal thatthe local unsteadl 
loads over most of the blade surface tend to suppress this motion, but those along Th“S 
surface from approximately 15% to approximately 45% of chord support the motion Z 
e ou -of-phase torsional motion, the pressure displacement function curves reveal that the 
local unsteady loads over the entire suction surface and most of the pr^ suit fendt 
suppress the in-phase torsional motion, especially those in the vicinity of the blade leading 

Similar results for bending vibrations are shown in Figure 9 The results nf th* f 

analyse for the in-phase bending vibration are in^d agr^IS “wh^h^ 
f the nonlinear Euler and the linearized potential analysis are in good agreement for the 
subresonant, out-of-phase motion. For the latter motion, the results of the two linear analyses 

m 7b"ade a cho1d er ' ?h eS t 0D ‘h' S “ Cti ° nS k urface from * blad ' '“<«■>* to approximately 

o!; the 1ft n^ of th I™ Z S nbt f ,ons m stabk . b «‘ local unsteady loads 
over the aft part of the suction surface are destabilizing for the in-phase motion 

Numerical results are given in Figure 10 indicating the behavior of the aerodynamic work 

nureT 6 Ver wi mter ^ n a t pl T e an § le for the subsonic NACA 5506 cascade undergoing 
pure torsiona (Figure 10a) and pure bending (Figure 10b) vibrations at unit frequency^ 

o?n tv, P f! 0113 ^ giV6n f ° r an interblade phase angle range from -90 dejto 

Phase" Lls e for r a ‘ IT ^ ^ 10 indicate the resonant interblade 

phase angles for a unit-frequency excitation. NPHASE, nonlinear, unsteady results are 

given for * = 0 deg, ±60 deg, ±90 deg, 120 deg, 180 deg, 240 deg a^d 270 deg; UN FLUX 

IZuhsh K S T ; int6rbIade phase an S le ’ an d a sufficient number of LINFLO 
results have been determined to essentially define W c as a continuous function of a The 

the th r analyses are in very g °° d a *~> 

at a - 120 deg, which is near the upstream resonance condition <x = <r+ = 117 12 de ff 

- ote that resonance phenomena are not accurately modeled by the current version of the 

SrjSfcir this anaIysi8 is appiied on a stretcM mesh with “« 

The work per cycle results in Figure 10 indicate that the single-degree-of-freedom tor- 
sional and bending blade motions at unit frequency are stable (W c < 0) and that the 

lower stability margin. They also LdLte leather cot 

d ffere i atur 7° f * lobal , aerodynamic response, which is associated with the 

different acoustic responses (i.e., subresonant, superresonant (1,0) etc ) that occur at dif 
feren. interblade phase angles. Note the abrupt changes in the gtoM tnite^y 

response behaviors, as predicted by the linearized analyses, at the interblade phase angles 
at which acoustic resonance occurs. P S 

th T genera1 ’ the resu ^ ts of the three analyses, LINFLUX, NPHASE and LINFLO for 
the subsonic unsteady flows considered in Figures 8 through 10, are in good agreement 
Differences exist, however, and the reasons for these are not clearly understood at the present 
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time. One possible cause is that 

composite H/C-mesh, respectively, are used i far . field Analytic, non-reflecting, Car- 

methods used to mod ^^ stea > he ^linearized calculations, whereas an axially stretched 
field conditions are applied dissinate outgoing disturbance waves. In future 

mesh is used in the nonlinear calculate t NPH ASE analysis so that 

^=^“t^r^ ta ‘^ iBd<riving ‘ he 

LINFLUX analysis from the nonlinear NP - 


Transonic Flow (M-oo - 0.8,0-00 58 deg) 


j ui At* cnrf arp Mach number distributions, ns 
The steady Mach number contours and W»^su^ MJ^ ^ Figure u _ with 

determined by the nonlinear Eu er ana ysis, CASPOF analysis. The Euler analysis 

surface Mach number distributions predicted by the CASPOF an y ^ ^ ^ ^ ^ 

predicts a mean flow exit Mach num er t +» upstream side and 

respectively. The Mach numbers at the base of of approximately 2.5% on the 
0.834 on the downstream side. Spurious P ent in the Eu i er predictions. The 

suction surface and 0.1% on the P« ss “' increased approximately 1.0% by the shock to 
total pressure loss on the suction * u ^ shock The po Jntial flow analysis predicts a mean 
approximately 3.5% downs rea respectively. The Mach numbers 

flow exit Mach number ^d ^ « “te s urf2 M^h numbers predicted by the 
at the base of the shock those on the suction surface just downstream 

two analyses are in good ag •> shock where small differences exist. 

of the leading edge and m the vicini y ’ hg transon i c NACA 5506 cascade are 

Unsteady aerodynamic respon! * P*e ing pure torsiona l and pure bending vibra- 

given in Figures 1- through out-of-phase component of the unsteady pressure 

tions at unit frequency, on ours o Tr^pLO and LINFLUX analyses, for an in-phase 
response, Im{p }, as predicte ^sing s hown in Figure 12. The pressure contours pre- 

(a = 0 deg) torsional blade vibration are ^ ™ * ^Lioa are in reasonably good 
dieted by the two analyses for this superreson (M) be noted that shock 

agreement, except in first . h armonic, unsteady 

fitting is used in the LINF y essures t'hat depend on shock displacement. The an- 

pressures and impulsive «n>» P p 12 shocks are captured in both the nonlinear 

“aSEiranalyJ so that tb. reSP °” SeS 

continuous, but undergo large changes in t e " c ‘“ y k per . cyc l e predictions for the tran- 
Pressure-displacement function ““tnTbl^e vTbratLs at w = 1 are shown 
sonic cascade undergoing in an ou ^ stable, and, according to the linearized 

in Figure 13. These superresonan a * indicated by the delta functions in Figures 13, 
potential analysis, the impulsive s shock impulse, predicted by LINFLO, 

contribute significantly to the stab, f m °" th ' e P sh ock and the shock 

displacement* 1 ^heXtcapt^ng uJ in the Euler calculations is based on the use of 
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directionaily-dependent flux evaluations. The directional differencing used in the nonlinear 
unsteady analysis is based on the nonlinear unsteady flow properties, whereas that used 
in the linear analysis is based upon the underlying steady flow properties. Although the 
time-accurate, nonlinear analysis predicts a pressure response with a sharply defined shock 
at each instant of time, shock effects are smeared in the pressure displacement function pre- 
dictions, because these are based on an integration in time of the local unsteady response 
over a complete cycle of blade, and therefore shock, motion. 

For the torsional motions considered in Figure 13, the results given by the nonlinear 
Euler analysis, the linearized Euler analysis, and the linearized potential analysis are in good 
agreement in subsonic regions, i.e., along the blade pressure surface and along the suction 
surface downstream of the shock. However, the results in the supersonic region upstream 
of the shock and in the vicinity of the shock show substantial differences, particularly those 
for the superresonant (1,0), out-of-phase, torsional blade vibration. Differences between 
the nonlinear and linearized Euler analyses are due in part to the different methods used 
to prevent disturbance reflections at the inflow and outflow boundaries, i.e, the use of a 
stretched mesh and the imposition of analytic far-field conditions, respectively, and to the 
different discrete approximations used in the vicinity of the shock. The differences between 
the results of the two linearized analyses are due primarily to the different methods (i.e., 
shock capturing and shock fitting) used in modeling unsteady shock phenomena. The shock 
capturing method currently used in the linearized Euler analysis is easy to implement, but 
it does not properly account for shock motion. Therefore, it does not lead to the prediction 
of a concentrated, anharmonic, pressure load at the shock. 

A similar set of results for bending vibrations is shown in Figure 14. These blade motions 
are also stable, even though the shock impulse, as predicted using LINFLO, is destabilizing 
for the out-of-phase motion. The results of the three analyses, particularly those of the 
nonlinear and linearized Euler analyses show significant differences in the vicinity of the 
shock and, for the in-phase motion, in the vicinity of the blade leading edge. The shock 
modeling effort in the linearized Euler analysis is continuing with the goals of achieving 
consistent and good agreement between the nonlinear and linearized Euler predictions in the 
vicinity of a shock. Also, an effort will be initiated to include analytic far-field conditions 
in the nonlinear code (NPHASE) so that the nonlinear and linear Euler calculations can be 
performed on the same mesh. This could eliminate the differences between the nonlinear 
and linearized Euler solutions at blade leading edges. 

Work per cycle versus interblade phase angle predictions are given in Figure 15 for the 
transonic NACA 5506 cascade undergoing pure torsional and pure bending blade vibrations 
at unit frequency. The three analyses, NPHASE, LINFLUX and LINFLO, provide results 
that show similar trends, and indicate that the unit frequency, bending and torsional vibra- 
tions of the transonic cascade are stable. 

7.3 Aerodynamic Excitation 

We proceed to analyze unsteady flows, excited by external aerodynamic disturbances, 
through the Tenth Standard Cascade operating at M = 0.5 and = 55 deg. The 

steady Mach number contours and surface Mach number distributions for the NACA 5506 
cascade, as determined by the NPHASE analysis, are shown in Figure 16 along with the 


65 



surface Mach number distributions and 39.9 deg, 

angle of 0.349 and 40.2 deg, respectively. 

Acoustic Excitations 

s b “,z. h — — * •> ■ ■ 

•riven bv (5.54) with d 2 m < 0, i.e., 

O ~ V 

Kj.^oo = i(=F|d T sol + M£oo&F«° COSn T“) ' 

The tangential wave number of the 

frequency, mterblade phase ang tra „ els aga inst the steady flow direction. 

S flowmonditions, greater mesh reso,ution is required to resolve 

quency w = 5 and mterblade phase “^000*4 to * waV elength of 3.73 blade chords, 
wave number magnitude is • . P lhe positive axial direction is 

and the propagation angle of Motions on the reference 

69.0 deg. The unsteady pressure field, and surtax :pr Utl anaJysis , are 

p^e e ;tXF^2 y ifL"he results of the two analyses are seen to be in excellent 

wavelength of 2.32 blade c , __ essure fields and pressure distributions on the 

positive axial direction is • S- _ na i V sis an d the linearized potential analysis are 

Made as predicted by the am again in very good 

presented in Figures 19 an • , c t r ~,m excitation is not as good as that 

U-nt, however wav^ength of the 

— - ■— “ ion errors to 

negligible levels. 

Entropic and Vortical Excitations 

Small-amplitude entropic and vortical disturbance s are convected at *££*££ 
steady background flow. A single angle, o, and the 

is decribed by specifiymg t e te ^ orl (2 ' 5) ' and (2 . 6 )]. For a vortical “citation 

r^U^tUple^tude of the normal component of the rotational velocity, 


i.e.. vr,_ oo • 


eN, where 


V/^—oo — i(#C_oo * C— co)/ 1^" 


(7.5) 
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directiotT SCtibi ? th ';° rtiCi ‘ y ' In' Tr‘ "*“** j ” «» M* free-stream 

direction, v* • e r , is determined by the divergence-free condition, i.e., -v R _ = Q 

e wave number k_oo, of a harmonic entropic or vortical excitation has a component 

"’7° 7 , , in the cascade tangential or ^-direction and a component k t = -ujV~ 1 = 

in the inlet free-stream direction. Therefore, _0 ° 

= ~(uj sec fi.oo -f aG~ l tan fi_ QO )e { + 

= -we r + (w tan + aG ~ 1 sec fi_ cc ) e ^ . 

To consider a gust that models the first harmonic component, i.e., the component at 
tlm Made passing frequency (BPF), of a wake excitation from upstream w. would set 

“ ,7 ’' 00 . Cz 7 ~27rG/G EXC and w = Vfexc = 2ir V exc /G E xc. Here G EXC is the 

blade spacing in the upstream row which moves at constant velocity V Exc e = V w e in 
the circumferential direction relative to the reference blade row, and 7-1 = -2 

“J tZT 7 TIi’T'T ° f ^ fundame " tal or first- harmonic disturbance. In 
the present study, an interblade phase angle of <r = -180 deg was chosen to represent the 

armomc component of a wake excitation from an upstream blade row with one-half the 
number of blades as in the reference blade row. The temporal frequency was set at „ = 5 
which corresponds to an upstream flow coefficient, i.e. ratio of mean flow axial speed to 
w ee speed, K-= cos fl^/W, of 0.3604. The wave number magnitude, |>c „| of this gust 

“-V s th “ 5 - 25 . “1 * ; direction of travel relative to the axial’ Llltion, 11 
a ~ tan ( / S,-oo/«£,-oo) ? is 36.5 deg. 

Results of the LINFLO and LINFLUX analyses for the convection of entropic [with 

(1 ’° d and vortlcaI Mh - ejv = (1,0)] gusts at BPF, i.e. u = 5 and <r - 

160 deg, through the NACA 5506 cascade operating at M_ TO = 0.5 and fL^ = 55 deg 

un^teT m , Fl§UreS l 1 thr ° Ugh 24 ' !t shouId be noted that the LINFLO predictions for the 
nsteady entropy and vorticity contours are based on closed form solutions in which these 

unsteady flow quantities are expressed as explicit functions of the drift and stream functions 
of the steady background flow (see [Ver93]). The predicted entropy fields for a prescribed 
entropic excitation at inlet are shown in Figure 21. The two linearized predictions for the 

lTnlrkedTlerT 111 It T ^ g °° d agreement ‘ However, the 

hTnumer^U show spurious distortions in the entropy contours, that are caused 

1 r Ca V°T’ m he , V,Cm,tieS ° f the blades and alon g blade wakes. A spurious 

Ztet° n NPHASE Py “l ?? 8teady background Predicted with the 

thro?g e h i "4 te Th t ,. i0n "“I 1 NACA 5506 CaSC4de is de P irted ” Figures 22 

in Tfr P«<fiffi vorticity and pressure fields, represented via contour plots of the 

he DmdimT°“T d VariaMeS ’ Me Sh ° Wn “ FigUres 22 and 23 ' respectively; and 

the predicted unsteady pressures acting on the surfaces of the reference blade, in Figure 04 

Mtte7n«rled Euler ' "T- ”7^' 22 7 “ g °° d 1 ualita «vo agreement, but the results 
of the Zest] w ^ 8 ° W SeVere “” S “ ‘he vorticity contours in the vicinity 

The ability to accurately predict the response of real blade cascades to entropic and 
vortical excitations using LINFLUX is currently impeded by the use of first-order aTu 
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rate, blade-surface, boundary conditions tott oApurious steady and 

unsteady Euler analyses. These Xfe surfaC es. Also, the entropy and 

first-harmonic, unsteady ent'opS . , ( m the steady background floiv are be 

;SSi^"i" » - vortical d.sturbances to further distort the 

“"^^^regoingli^ations,^ — 

LINFLUX analysis differ from those predated Lby W a strong impact on the 

the predicted vorticit, feMs <**“ 23) “‘f 

resulting unsteady pressure respo linearized Euler and potential analyses or 

pressure distributions (Figure 24), predict y ^ ^ ^ agreement Nevert heless. the 

the vortical gust at w 0 an a unsteady flows excited by vortical disturbances must 
issues should be addressed m future wor . 
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8. Concluding Remarks 


to rnh and f fl ? Cient ’ three ' dlmens ional, unsteady, fluid-dynamic analyses are needed 
to enhance our abihty to understand and predict the aeroelastic (flutter and forced vt 

ion) and aeroacoustm (noise generation, transmission and reflection) characteristics of 

^ dmg - SuC i analyses must account for nonisentropic and rotational invis 
d effects to model swirling flows with strong shocks and, eventually, for viscous effects to 

model separated flows Motivated by these needs, independent, but coordinated ^search 
programs are being conducted at NASA Lewis Research Center Mississinoi St^e TT f 

rusted 6 Unit 1 TeChn0l0g T RCSearCh Center t0 d -lo P a^ 

the presLTuTRcTIo^ na ^ f ° r ^ r ° elast * c and a«roacoustic design applications. Under 
. P resent UTRC program, work has been focused on formulating linearized invisrid „ n - 
s eady aerodynamic and numerical models, based on the Euler equations of fluid m f ’ f 

moTelf h men f nal thr ° U ? Sin§le bkde r ° WS ° PeratiDg witbin ^rical duTs These 

rad^LIN^LU^^nd^the'LI^FLUX^^'^ ^h tlenS *° n ^' '~ d ' ^ 

unsteady ctcadelwl ' ^ ^ ”*** ™ ^ benchmark 

In formulating the linearized unsteady aerodynamic model, the unsteady flow is regarded 
as a small perturbation of a nonlinear background flow that is steady fn the blade row 
frame of reference. This leads to time-independent, nonlinear, governing equations for the 
steady background flow, and time-independent, linear, variable-coefficient, equations for the 
first-harmonic unsteady perturbation. The variable coefficients in the unsteady equation! 
of ^ background flow - This of linearization accounts for t he effects 

on the 4 f 6 g r metr ^ mean blade loading ’ mean swirl str ong shocks and their motions 
exdtations y yMm,C resp ° nses to P^ribed structural and external aerodynamic 

of , th ' 1 ■ T Dt f ? rmulat !°“; both the s ‘eady and unsteady equations are expressed in terms 

of spat, al independent vanables that describe the reference or steady-state LsTuons of the 

are sT A P °* ntS that make U P the P h y sica l solution domain. Thus, both sets of equations 
e solved on a fixed grid in reference physical space to determine the values of the steady and 
unsteady flow variables at the mean and instantaneous positions, respectiv^ of the ^Wng 

nhv S P °iH S TtlS appr ° ach f HC93a > Gl193 ) in which the time-dependent deformation of the 
physical domain is represented by source terms in the unsteady field equations simnlifie* th 

oin rt l i bound « y condi«Lns ld rsr;s^ £ 

account for the effects of mean swirl and the radial dependence of the steady and unsteady 

steld a ? Stl be developed to cor nplete the formulation of the three-dimensional 
steady and first-harmonic, unsteady flow problems. aimensional 

, T he fl nUmerkal mod els developed herein for determining the steady and linearized un- 
steady flows are based on the nonlinear, time-accurate, Euler analysis developed in [WJS881 

"Z^rT T dent tranSf0rmati0n fr ° m reference " comp u raS CO-' 

dinates and a pseudo-time marching procedure in which the steady and the first-harmonic 
unsteady °w SS arc determined by marching implicitly in pseudo-time, using S ttaTstep 

approximated S ° IUti ° nS "* delermined - The Psuudo-time derivatives are 

pproximated by first-order difference approximations, and, except for this chance the 

men cal model for the steady flow is essentially identical to the one used in [WJSSSJ for 
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nonlinear unsteady flows. The linearized unsteady equations contain source terms, that 
depend on the blade motion, and the spatial differencing directions used to evaluate the 

“Tpr y evto^ m ° del ° f lW f 8 V "* 

JTTLTIU multi-block, finite- volume c ofi®i caUed NPHASE to p- , a-- 
dimensional unsteady aerodynamic analysis for vibrating cascades The NPHASE code has 
been used al the bas^s for implementing the steady and linearized unsteady aerodynamic 
models described in this report. In particular, nonlinear steady flows are predicted using a 
slightly updated version of NPHASE, and linearized unsteady flows are determmed usmg a 
new code called LINFLUX, which has been constructed, from the original NPHASE anal- 

7Z Aspartofthe present effort, the surface-boundary-condition implementation used in 

NPHASE has been modified to reduce spurious numerical losses, and one- and approximate 
two-dimensional, far-field conditions [Gil90] have been incorporated 

unsteady flows excited by external aerodynamic disturbances. Unfortunately, the latt 
conditions have been found to result in significant reflections, at the inflow and outflo 

boundaries, of outgoing unsteady disturbance waves. Th^^Uons 

far-field conditions have been constructed and implemented into LINFLUX. lhese cond 
seemto eliminate the non-physical reflections, indicating that similar conditions should be 

T " - -e appu«l u to 

transonic flows through the Tenth ^ ^ acoustic disturbim ces at inlet 

tTZ 3” inlet. For each example, LINFLUX pre- 

dktions have been compared with other analytical predictions; namely predictions based 

nonlinear, analysis NPHASE, and/or those based on the linear, zed 

POt Th“st LINFLUX gives response information that is in 

The results ot t y earlier analyses NPHASE and LINFLO, for unsteady 

close agreement with dWurbances that travel towards the 

subsonic flows exci y INFLUX capabilities are needed for transonic flows and for 

blade toW t Howeve, ,mproved Modeling or directional differencing, 

flOW d in LINFLUX should be modified so that LINFLUX and NPHASE will provide consis- 
used m LINFLUX should bem ^ for unsteady transonic flows. Also the 

transonic and vortical gust calculations. 
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The primary reason for developing nonlinear and linearized unsteady Euler analyses for 
blade rows is the need to predict three-dimensional unsteady flows in which the effects of 
mean swirl and radial variation are important. Thus, the follow-on to the present research 
effort will be directed towards implementing the aerodynamic and numerical formulations 
lined in this report, into a three-dimensional version of the LINFLUX code. In addition 
to the problems associated with the numerical modeling of shocks and vortical gusts, another 
major difficulty must be overcome. This involves the development and implementation of 
appropriate far-field conditions for three-dimensional steady and unsteady flow calculations. 
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Figure 1: Rotating axial compressor blade row operating within a cylindrical duct. 
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Figure 2: Two-dimensional section of a ““^oTthetuctil luXt of the nth 

& ” "r t - - — - of tte nth 








(a) 



(c) 



F, gate 4: Computational 

state grid in physical space; and (c) un.torm 
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Figure 5: H-meshes used for steady (NPHASE) and linearized unsteady (LINFLUX) Euler 
calculations: (a) subsonic flow; (b) transonic flow. 
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Figure 7: Contours of the out-of-phase component of the unsteady pressure, i.e., Im{p}, for 
the subsonic (A/_ <*, = 0.7, = 55 deg) NACA 5506 cascade undergoing an out-of-phase 

torsional blade vibration about midchord (a = 2 deg, u = 1, a = 180 deg): (a) linearized 
potential analysis (LINFLO); (b) linearized Euler analysis (LINFLUX). 
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r^i^rchola at':= /d^d „-!:(.) **- (, - 0 deg) — biade 
motions; (b) out-of-phase (<r = 180 deg) motions. 
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Figure 9. Pressure-displacement function distributions and work per cycle predictions for the 
subsonic (Af_oo = 0.7, f!_oo = 55 deg) NACA 5506 cascade undergoing bending vibrations 
at h y = 0.01 and u = 1: (a) in-phase (a = 0 deg) bending motions; (b) out-of-phase (a = 
180 deg) motions. 
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Figure 10: Work per cycle vs. 

55 deg) NACA 5506 «»«*• -^^j£2r,£V = 0 01. 


55 deg) NAUA ooUD cascauc un^ 6 — & , , __ n nl 

aLut midchord with a = 2 deg; (b) bending v.bratrons w,th h, - 0.01. 
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Figure 11. Mach number contours and blade surface Mach number distributions for steady 

transonic flow at = 0.8, and = 58 deg through the NACA 5506 (10th Standard) 
cascade. 
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analysis (LINFLO); (b) linearized Euler analysis (LIN )• 
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Figure 13: Pressure-displacement function distributions and work per cycle predictions for 
the transonic (M_ = 0.8, fi_ = 58 deg) NACA 5506 cascade undergoing torsional blade 
vibrations about midchord at a = 2 deg and u> = 1: (a) in-phase (a = 0 deg) torsional blade 
motion; (b) out-of-phase (<r = 180 deg) motion. 


91 






Figure 14: Fr^re-d^placement funcUon “tXuudergTug Ending vibra- 

^“and Stt W in-phase «r = 0 deg) bending vibration; (b) out-of-phase 
( a = 180 deg) bending vibration. 
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Figure 15: Work per cycle vs. interblade phase angle for the transonic (M- x = 0.8 fl_ = 
58 deg) N AC A 5506 cascade undergoing blade vibrations at u, = 1: (a) torsional vibrations 
about midchord with a = 2 deg; (b) bending vibrations with h y = 0.01. 
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Figure 17: Contours of the in-phase component of the unsteady pressure due to the interac- 
tion of an acoustic excitation from upstream (?/._«, = (1,0), w = 5, <7 = -90 deg) with the 

r bS r T C Ti M r r °' 5, = 55 deg) NACA 5506 cascade: ( a ) LINFLO calculation; (b) 

LUNbLUX calculation. 
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Figure 18: Unsteady surface pressures tnbuti. (M-T 

excitation from upstream (p/,-oo l ’ h 
n 5 fl = 55 deg) NACA 5506 cascade. 






Figure 19: Contours of the in-phase component of the unsteady pressure due to the interac- 
tion of an acoustic excitation from downstream (p/, +00 = (1,0), u = 5, a = -90 deg) with 

^t S i?T n nv (M rr 0 - 5 ’ ^ = 55 deg) NACA 5506 CaSCade: < a ) LINFL0 calculation; 
(bj UINr LUX calculation. 


97 



Figure 20: Unsteady surface pressure di 




Figure 21: Contours of the in-phase component of the unsteady entropy for an entropic gust 
(s_oo = (1,0), uf = 5, a = — x) interacting with the subsonic (M_ <*> — 0.5, f!_oo = 55 deg) 
NACA 5506 cascade: (a) LINFLO calculation; (b) LINFLUX calculation. 
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Figure 22: Contours of the in-phase component of the unsteady vorticity for a vortical gust 
. enr = (1 01 u = 5, cr = -tt) interacting with the subsonic (Af-co - - 00 

deg) NACA 5506 cascade: (a) LINFLO calculation; (b) LINFLUX calculation. 
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Figure 23. Contours of the in-phase component of the unsteady pressure response for a 
vortical gust = (1, 0), u = 5, a = — 7r) interacting with a subsonic (M_ <*, = 0.5, 

oo ^ deg) NACA 5506 cascade: (a) LINFLO calculation; (b) LINFLUX calculation. 
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Figure 24: Unsteady surface pressure 

(VB..»-e w = (l,0),«=*5, <r = ~*) 

5506 cascade. 


distributions due to the interaction 
and the subsonic (M_ «> = 0.5, ft-co 


of a vortical gust 
= 55 deg) NACA 
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